PCX 



world iothxictual property organization 

international Bureau 




imgRNA-noNAL AP^T^jnj^^^^ PATENT co 0PERATION TOEATY 

International Patent Ckuxifi™.tY». 6 . * I — __ 

WO 99/26078 

27 May 1999 (27.05.99) 




(51) International Patent Classification 6 
COIR 33/24, 33/16 



AI 



(J 1) International Publication Number; 
(43) International Publication Date: 



(21) InternatioDal Application Number: PCT/US98/24220 
| (22) International Pfling Date: 13 November !998 (13.1 1.98) 



J (30) Priority Data: 
6G7065.318 
Not furnished 



13 November 1997 (I3.ti.97) US 
12 November 1998 (12.11.98) US 



(71X72) Applicant and Inventor: MILLS, RandeM, U (US/USl* 
I 1780 Greenbriar Court, Yardlcy, PA 19067 (US). 

I U4) Ag ?^ ^9°™' Gil *°. A. Jr. ct ah; Lahive A Cockfiefd 
I LLP, 28 Stale Street, Boston, MA 02109 (US). 



(81) Designated Slates: AL, AM, AT, AU. AZ> BA BB BG fiR 

nr. ca, ch. cn, cu, cz. de. dk?3^ es k cb cd 

KZ, LC, LK, LR, LS, LT. LU, LV. MD, MG MK MN 
MW. MX, NO, NZ, P U PT. RO. RU, SD SE . kc sia 
SL. TJ, TM, TR, TT. l/A. UG, UZ. VN. W ^ I 
patent (GH CM, KE, LS, MW, SD. SZ> UG, ZW), Eurasian 
patent (AM, AZ, DY, KG. KZ, MD, RU, TJ r TM) European 
patent (AT. BE. CH, CY, Da DSC, ES, FI, FR, GBGR 
IE. IT. LV. MC. NL, PT, SE), OAP1 patent (BP BJ CF 
CO, CI. CM. GA. GN. GW, ML, MR. NE, SN, TD, TG) 

Published 

With international search report. 



COPY 



(S4) Ttllc: RESONANT MAGNETIC SUSOiPTIBILITY FMAGING (RcMSI) 
(57) Abslrsct 



f/Sa - 



HcoMsrsocnoM 

WCiSSOi 



n 






UCMTOi 




- ~U6b 



A three-dimensional 
I magnetic susceptibility map of an 
object including a patient placed 
in a magnetic field is generated 
from a tKrec-dimcnstonal map 
I of the resonant radio frequency 
J (RF) magnetic flux external to the 
j patient. The magnetic susceptibility 
Of each voxel is (fclcrmined from 
I the shift of the Larmor frequency 
J due to the presence of the voxel in 
the magnetizing field. The intensity 
I variation of the transverse RF field 
J over space is used to dcfxnninc 
the coordinate location of each 
voxel. The RF field is the near 
field whichr is a dtpole that serves 
as a basis dement to form a unkjue 
reconstruction. The geometric 
system function corresponding to 
a dipofe which determines the 
spatial Intensity variations of the 
RF field is a band-pass for kp~k x . 
|n the limit, each volume element 
is reconstructed independently 
in parallel with all other volume 
elements such that the scan time 
[ is *> Sweater than the nuclear free 

iS^^ *. La™ n m) fmm m ^ 

pm* navel magnetic kugtag (MRI) invention X B o^l7 ^ 'J**^ Wi * ^ uni 1 oe , 

vomica. imae „ an, imag « _ £££ SntS^: i^"'^ ' 



FOR THE PURPOSES OF INFORMATION ONLY 



Codes used lo identify States party to the PCT on the front pages of pamphlets publishing international applications under the PCX. 



AL 




RS 


Spun 


LS 


Lesotho 


St 


Slovenia 


AM 


Anncrna 


Ft 




LT 


Lithuania 


SK 


Slovakia 


AT 


Austria 


F« 


France 


UJ 


Luacmbocrrg 


SN 


Senegal 


AU 


Australia 


CA 


Gabon 


LV 


Latvia 


SZ 


Swaututd 


KZ 


A&citxaijan 


CB 


United Kiogdom 


MC 


Monaco 


Tf> 


Chad 


BA 


Bosnia and lferzegovttit 


CE 


Georgia 


MD 


Republic of Moldova 


TG 


Togo 


BB 




CH 


Ch«na 


MC 


Madagascar 


TJ 


Tajikistan 


B6 


Belgium 


CN 


Guinea 


MK 


1 Tie forma Yugoslav 


TM 


Turkmen tsl an 


BF 


Durkraa Faso 


CK 


Greece 




Republic of Macedonia 


TR 


Turkey 


BG 


Bulgaria 


UU 


Hungary 


Ml, 


Mali 


TT 


Trinidad and Tobago 


BJ 


Benin 


IE 


Ireland 


MN 


Mongolia 


UA 


UklBroc 


BR 


Dnza 


It, 


Unci 


MR 


Mauritania 


UC 


Uganda 


By 


BcttfM 


IS 


Iceland 


MW 


Malawi 


us 


Untied Stales of Am r 


CA 


Canada 


IT 


Italy 


MX 


MeaifO 


uz 


Ucbekbtaa 


o 


Centra) African Republic 


JP 


Jap in 


NK 


Nigw 


VN 


Vid Nam 


cc 


Congo 


KE 


Kenya 


NL 




TU 


Yugoslavia 


ai 


Switzerland 


KG 


Kyrgyzstan 


NO 


Norway 


zw 


Zimbabwe 


CI 


Oe dl*ci«c 


KP 


Democratic People" i 


NZ 


New Zealand 






CM 


Cencrooo 




Republic of Korea 


PL 


Poland 






CM 


Cbint 


KR 


Republic of Korea 


rr 


Portugal 






CU 


Cuba 


KZ 


KiLaksXan 


RO 


Rocntfiia 






CZ 


Oech Rept&Kc 


LC 


Skim Lucia 


RU 


' Russian Federation 






Dt 


Gmni/iy 


U 


LicchtcMiehi 


st> 


Sudan 






VK 


DenmaA 


LK 


Sri Lenta 


SE 


Swede* 






EK 


Estonia 


LR 


Lfce/h 


$C 


Singapore 







WO 99/26078 PCT/VS98/24220 



10 



20 



RESONANT MAGNETIC SUSCEPTIBILITY IMAGING (ReMSI) 

BACKGROUND 

This invention relates to method and apparatus for imaging a 
body. In particular, the invention provides a magnetic susceptibility 
image of an animate or inanimate body. Many imaging techniques 
exploit some natural phenomenon which varies from tissue to tissue, 
such as acoustic impedance, nuclear magnetic relaxation, or x-ray 
attenuation to provide a contrast image of the tissue. Alternatively, 
some imaging techniques add a substance such as a positron or 
gamma ray emitter to the body and construct an image of the body 
by determining the distribution of the substance. Each imaging 
technique possesses characteristics which result in certain advantages 
15 relative lo other imaging Ie chniqucs. For example, the short imaging 
nmc of x-ray contrast angiography reduces motion artifacts. In 
addition, ihe high resolution of x-ray contrast angiography renders 
this technique superior to many prior known imaging techniques for 
high resolution imaging of veins and arteries. However, x-ray 
contrast angiography is invasive, requires injection of a noxious 
contrast agent, and results in exposure to ionizing radiation. Thus, it 
is not typically employed except for patients with severe arterial or 
venous pathology. 

25 Nuclear Magnetic Imaging (NMR) which is commonly called 

magnetic resonance imaging (MRI)) emails magnetizing a transverse 
tissue slice with a constant primary magnetic field in a" direction 
perpendicular to the slice, and further magnetizing the slice by 
applying a gradient in the plane of the slice. A radiofrequency pulse 

30 excites selected nuclei of the slice. The excited nuclei relax and emit 
energy, i.e.. radio signals, at frequencies corresponding to local 
magnetic fields determined by the gradient. A Fourier analysis of the 
emitted signals provides the signal intensity at each frequency, 
thereby providing spatial information in one dimension. Repeating 

3:> the excitation of the nuclei and obtaining the Fourier spectrum of the 
emitted signals, as the gradient rotates in the plane of the slice, 
provides a two-dimensional image. 

MRI is of primary utility in assessing brain anatomy and 
40 pathology. But long NMR relaxation times, a parameter based on how 
rapidly excited nuclei relax, have prevented NMR from being of utility 
as a high resolution body imager. The most severe limitation of NMR 
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technology is that for spin echo imaging n, the number of free 
induction decays ("FIDs"), a nuclear radio frequency energy emitting 
process, must equal the number of lines in the image. A single FID 
occurs over approximately 0.1 seconds. Not considering the 
5 spin/lattice relaxation time, the time for the nuclei to reestablish 

equilibrium following an RF pulse, which may be seconds, requires an 
irreducible imaging lime of n times 0.1 seconds, which for 512 x 512 
resolution requires approximately one minute per each two 
dimensional slice. This represents a multiple of 1500 times longer 

10 that the lime that would freeze organ movements and avoid image 
deterioration by motion artifact. For example, to avoid deterioration 
of cardiac images, the imaging time must not exceed 30 msec. A 
method for speeding NMR imaging flips the magnetization vector of 
the nuclei by less than 90 degrees onto the x-y plane, and records less 

15 FIDs. Such a method, known as the flash method, can obtain a 128 x 
128 resolution in approximately 40 seconds. Another technique used 
to decrease imaging time is to use a field gradient and dynamic phase 
dispersion. corresponding to rotation of the field gradient, during a 
single FID to produce imaging times typically of 50 msec. Both 

20 methods produce a decreased signaUo-noise ratio ("SNR H ) relative to 
spin echo methods. The magnitude of the magnetization vector which 
links the coil is less for the flash case because the vector is flipped 
only a few degrees into the xy-plane. The echo-planar technique 
requires shorter recording times with a concomitant increase in 

25 bandwidth and noise. Both methods compensate for decreased SNR 
by increasing the voxel size with a concomitant decrease in image 
quality. Physical limitations of these techniques render obtaining 
high resolution, high contrast vascular images impractical. 

30 Thus, it is an object of the invention to provide high resolution 

multi-dimensional images of tissue. 

It is another object of the invention to provide multi- 
dimensional magnetic susceptibility images of an object. 

It is yet another object of the invention to provide high 
resolution multi-dimensional images of the cardiopulmonary system. 

It is yet another object of the invention to provide a magnetic 
40 susceptibility image of a body. 
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SUMMARY OF THE INVENTION 
These and other objects of the invention are attained by providing an 
5 apparatus for obtaining a multi-dimensional susceptibility image of a body. Hie 
apparatus includes a radiation source for magnetizing the body with a magnetic 
component of a first radiation field. The apparatus also includes a first detector for 
measuring the magnetic component of the first radiation field in the absence of the 
body in a volume to be occupied by the body. The apparatus further includes a 

10 source for applying a second radiation field to the body, to elicit a third radiation 
field from the body. A second detector senses this third radiation field, and produces 
a signal that a reconstruction processor employs to create the magnetic susceptibility 
image of the body. 

One practice of the invention provides a method for determining the 

IS distribution of radiation within a magnetized body, emanating from the body in 
response to an excitation radiation. The method includes the steps of measuring the 
emanated radiation over a three-dimensional volume by an array of detectors, and 
uniquely correlating each frequency component of the detected radiation with 
locations within the body producing that frequency component. 

20 The invention is in part based on the realization that matter having a 

permeability different from that of free space distorts a magnetic flux applied 
thereto. This property is called magnetic susceptibility. An object, herein called a 
phantom, can be considered as a collection of small volume elements, herein 
referred to as voxels. When a magnetic field is applied to the phantom, each voxel 

25 generates a secondary magnetic field at the position of the voxel as well as external 
to the phantom. The strength of the secondary magnetic field varies according to the 
strength of the applied field, the magnetic susceptibility of the material within the 
voxel, and the distance of the external location relative to the voxel. For example, 
U.S. Patent No. 5,073,858 of Mills, herein incorporated by reference including the 

30 references therein, teaches that the net magnetic flux at a point extrinsic to a phantom 
to which a magnetic field is applied, is a sum of the applied field and the external 
contributions from each of the voxels. The f 858 patent further teaches sampling the 
external flux point by point and employing a reconstruction algorithm, to obtain the 
magnetic susceptibility of each voxel from the sampled external flux. 

35 Unlike the '858 patent that relies on a static response from a 

magnetized body to determine the magnetic susceptibility of the body, a preferred 
practice of the invention elicits a radiative response from a magnetized body by 
subjecting the body to a resonant radiation field. One embodiment of present 
invention generates three-dimensional magnetic susceptibility image of an object 

40 including a patient placed in a magnetic filed from a three- 
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dimensional map of a radio frequency (RF) magnetic field external to 
the patient, induced by subjecting selected nuclei of the body to a 
resonant RF field. Application of an RF pulse to the body causes the 
body to emit the RF magnetic flux external to the body. A Fourier 
5 transform of this external flux produces its frequency components 
("Larmor frequencies"). Each Larmor frequency determines the 
magnetic susceptibility of the voxels of the body producing that 
Larmor frequency. Further, the intensity variation of the external RF 
field over a three-dimensional volume of space is used to determine 
10 the coordinate location of each voxel. 

The radiation source for magnetizing a body to be imaged can 
be a direct current ("DC") magnet, including a superconducting 
magnet. The radiation sources and amplifiers for applying an RF 

15 pulse to the magnetized body are well known in the art, and include, 
but are not limited to 5 klystrons, backward wave oscillators, Gunn 
diodes, Traveling Wave Tube amplifiers. A preferred embodiment of 
the invention employs superconducting quantum interference devices 
(SQUID) as detectors for sensing the external RF field. SQUIDs 

20 advantageously allow nulling the applied magnetic flux in order to 

measure small external fields produced by precessing nuclei. Because 
the contributions of voxels to an external field at a detector typically 
drops off as the cube of the distance between the voxel and the 
detector, some embodiments of the invention employ a magnetizing 

25 field whose amplitude increases as the cube of the distance from the 
detectors. Such a spatial variation of the magnetizing field "levels" 
the magnitude of the external RF radiation from voxels at different 
locations relative to the detectors. 

30 One practice of the invention obtains a . three-dimensional 

magnetic susceptibility map of a magnetized body from a three- 
dimensional map of a secondary magnetic f]ux produced by the 
magnetized body, and detected over a three-dimensional volume of 
space external to the body, herein referred to as the "sample space." 

35 The extrinsic magnetic flux is sampled at least at the Nyquist rate, i.e., 
at twice the spatial frequency of the highest frequency of the Fourier 
transform of the magnetic susceptibility map of the phantom, to allow 
adequate sampling of spatial variations of the external magnetic flux. 
This practice of the invention preferably employs a Fourier transform 

40 algorithm, described in Fourier Transform Reconstruction Algorithm 
Section, to form the magnetic susceptibility map of the object" 
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One embodiment of the present invention emplovs nuclear 
magnetic resonance (NMR) to induce a magnetized body to emit an 
external radiation having a magnetic field component. In particular, 
application of an RF pulse, resonant whh selected nuclei of a 
5 magnetized body, can polarize the nuclei through rotation of their 

magnetic moments. The polarized nuclei within a voxel precess about 
the local magnetic field in the voxel at a Larmor frequency 
determined in part by the applied magnetic field at position of the 
voxel and the susceptibility of the voxel. The superposition of 
10 external RF fields produced by all the voxels of the body creates the 
total external RF field. Thus, the external RF field contains frequency 
components corresponding 10 precession frequencies of the nuclei in 
different voxels of the body. 

1 5 One practice of the invention detects the external RF field in the 

near field region where the distance of a detector sensing radiation 
from a voxel at a distance r from the detector is much smaller than 
the wavelength A of the radiation emitted by the voxel, i.e.. r«X (or 
kr«\). The near fields are quasi-stationary, that is they oscillate 

20 harmonically as c"", but are otherwise static in character. Thus, the 
transverse RF magnetic field of each voxel is that of a dipole. In one 
embodiment, an array of miniature RF antennas sample the external 
RF field over a three-dimensional volume of space that can be either 
above or below the object to be imaged. The frequency components 

-5 or the detected RF signals determine the magnetic susceptibility of 

voxels that give rise to the RF signals, and the location of each voxel is 
determined through measurements of the spatial variations of the 
intensity of the external RF field at a given frequency. Thus, the 
frequency components of the external RF radiation, and the intensity 

30 vanations of the external RF radiation provide the necessary 
information for providing a magnetic susccptibilitv map of the 
magnetized phantom, such as a human body. Such a susceptibility 
map can be employed to obtain anatomical images of a human body 

^ based on selected physiological parameters. 

Another practice of the invention employs paramagnetic and/or 
diamagnetic substances present in a body to be imaged to cause 
variations of local magnetic fields in different parts of the body, 
thereby providing a susceptibility image of the body. In particular 
40 an excitation field, such a magnetic field, polarizes the paramagnetic 
and/or diamagnetic substances of the body. The polarized substances 
contribute to the local magnetic field at the positions of the voxels 
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comprising the body. The amount of this contribution varies from one 
voxel to another depending on the variation of the concentration of 
the substances throughout the body. An excitation of selected nuclei 
of the voxels cause the nuclei to provide an external RF field through 
5 precession about the local magnetic fields. The external RF field 

includes a range of frequencies relating to the local magnetic field at 
position of each voxel. Thus, an analysis of the frequencies lead to 
information regarding the distribution of the paramagnetic or 
diamagnetic substances throughout the body. 

10 

Magnetic Susceptibility 

All matter has a permeability different than that of free space 
thai distorts an applied magnetic flux. This property is called 

15 magnetic susceptibility. An object, herein called a phantom, can be 
considered as a collection of small volume elements (hereafter called 
"voxels"). When a magnetic flux is applied to the phantom, each voxel 
generates a secondary magnetic flux at the position of the voxel as 
well as external to the phantom. The strength of the secondary 

20 magnetic flux varies according to the strength of the applied flux, the 
magnetic susceptibility of the material within the voxel and the 
distance of the external location relative to the voxel. In one 
embodiment described in Mills (I] which is herein incorporated by 
reference including the references given therein in the Reference 

25 Section, the net magnetic flux at a point extrinsic to the phantom is a 
sum of the applied flux and the contributions of each of the voxels 
within the object. This flux is point sampled over a three dimensional 
space, and the magnetic susceptibility of each voxel is solved by a 
reconstruction algorithm. 

30 In an embodiment of the present invention, the three- 

dimensional magnetic susceptibility map of an object including a 
patient placed in a magnetic field is generated from a three- 
dimensional map of the transverse resonant radio frequency (RF) 
magnetic fiux external to the patient. The magnetic susceptibility of 

35 each voxel is determined from the shift of the Larmor frequency due 
to the presence of the voxel in the magnetizing field. The intensity 
variation of the transverse RF field over space is used to determine 
the coordinate location of each voxel. The RF field is the near field 
which is a dipole that serves as a basis element to form a unique 

40 reconstruction. The geometric system function corresponding to a 
dipole which determines the spatial intensity variations of the RF 
field is a band-pass for In the limit, each volume element is 
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(2) 



reconstructed independently in parallel with ail other volume 
elements such that the scan time is no greater than the nuclear free 
induction decay (FID) time. 

5 Secondary Magnetic Field 

The magnetic moment, of each voxel within the phantom is 
given by the product of its volume magnetic susceptibility, ^, the 
magnetizing fiux oriented along the z-axis. B , and the volume of the 
10 voxel, AV. 

The magnetic moment of each voxel is a magnetic dipole. And the 
phantom can be considered to be a three-dimensional array of 
magnetic dipoles. At any point extrinsic to the phantom, the z- 
15 component of the secondary flux. B. from any single voxel is 

B =//!_ ; ~ 

where x, y, and z are the distances from the center of the voxel to the 
sampling point. It is shown in APPENDICES MV that no geometric 
distribution of magnetic dipoles can give rise to Eq, (2). Therefore, 

20 the flux of each magnetic dipole (voxel contribution) forms a basis set 
for the flux of the array of dipoles which comprise the magnetic 
susceptibility map of the phantom. 

Eq. (2) is a system function which gives the magnetic flux 
output in response to a magnetic dipole input at the origin. The 

25 phantom is an array of spatially advanced and delayed^ dipoles 

weighted according to the magnetic susceptibility of each voxel; this is 
the input function. The secondary flux is the superposition of 
spatially advanced and delayed flux, according to Eq. (2); this is the 
output function. Thus, the response of space to a magnetized 

30 phantom is given by the convolution of Eq. (2) with the series of 
weighted, spatially advanced and delayed dipoles representing the 
susceptibility map of the phantom. 

In Fourier space, the output function is the product of the 
Fourier transform (FT) of the system function and the FT of the input 

35 function. Thus, the system function niters the input function. The 

output function is the flux over all space. However, virtually all of the 
spectrum (information needed to reconstruct the magnetic 
susceptibility map) of the phantom exists in the space outside of the 
phantom because the system function is essentially a band-pass filler. 

40 This can be appreciated by considering the FT. //[*,.*.], of Eq. (2): 
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where k p is the spatial frequency in the xy-plane or Vplanc and k - 'S 
the spatial frequency along the z-axis. H{k,X\ is a constant for k and 
*. essentially equal as demonstrated graphically in FIGURE 1c. 



Band-Pass Filter 

A magnetic field with lines in Ihe direction of the z-axis applied 
to an object comprised of magnetically susceptible material 
magnetizes the material so that a secondary field superposes the 
applied field as shown in FIGURE 9. The secondary field outside of 
the object (phantom) and detected at a detector 301 is that of a series 
of magnetic dipolcs centered on volume elements 302 of the 
magnetized material. In Cartesian coordinates, the secondary 
15 magnetic flux, B', at the point (x, y ; z ) due to a magnetic dipole having a 
magnetic dipole moment m. at the position (x^.v,.^) is 

B (5) 
20 where i. is the unit vector along the z-axis. The field is the 

convolution of the system function, h(x.y.z) or h{p.<p. z ) (the left-handed 
part of Eq. (5)), with the delta function (the right-hand pan of Eq. 
(5)), at the position U OT y 0 ,Co). A very important theorem of Fourier 
analysis states that the Fourier transform of a convolution is the 
25 product of the individual Fourier transforms [2). The Fourier 
transform of the system function, h(x t v,z)orh(p,<t> xZ ) t is given in 
APPENDIX V. 

The z-component of a magnetic dipole oriented in the z- 
direction has the system function, Afx.v.j). which has the Fourier 
30 transform. which " shown in FIGURE lc. 
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The output function, the secondary maenetic field is the 
convolution of the system function, h(x.r.z)~iht geometric transfer 
function for the z-componcnt of a z-oricnted magnetic dipole with the 
input funciion-a periodic array of delta functions each at the position 
of a magnetic dipolc corresponding to a magnetized volume element 

(2z--.r-y) - 

i^J + v 2 + : 3 ) 5 '-' 2L m : S (*- ftx...v-nv , z -n z .) (8 ) 
The Fourier transform of a periodic array of delta functions (the 
nght-hand side of Eq. (8» is also a periodic array of delta functions in 

K - S p 3 C C . 



By the Fourier Theorem, the Fourier transform of the spatial output 
function, Eq. (8), is the product of the Fourier transform of the system 
function g,vcn by Eq. (7). and the Fourier transform of the input 
function given by Eq. (9) 



v 



In the special case that 

, _ *' = * (11) 

the Fourier transform of the system function (the left-hand side of Eq 
(10)) is given by H ' 

Thus, the Founer transform of the system function band-passes the 
Founcr transform of the input function. Both ,hc input function (the 
nght-hand part of Eq. (8)) and its Fourier transform ((he ri e hl-hand 
part of Eq. (10)) are a periodic array of delta functions No " 
frequencies of the Fourier transform of the input function are 
attenuated: thus, no information is lost in the case where Eq (11) 
holds. Thus, the resolution of ihe reconstructed magnetic 
susceptibihty map is limited by the spatial sampling rate of the 
secondary magnetic field according ,o the Nyquist Sampling Theorem. 

Reconstruction of the Magnetic Susceptibility Map 

Superconducting Quantum Interference Devices (SQUIDs) allow 
the user to null the ambient or applied magnetic flux and to measure 
very small secondary magnetic fluxes ( lO"" Testa) The flux 

tvmcaHv'm-T T' ° f hUmaa ,iSSUe of mi,Hmctcr dimensions is 
«yp.cally 10 Tesla. For example, representative parameters of Eq. (1) 
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are * = !0^, B_ = 10 s G % and AV = loW which results in m, = 10"* Gcm\ 
Substitution of this voxel magnetic moment into Eq. (2) with 
x = 0;y = 0;* = 10c/n results in B = 10~ 5 G = \ nT . Eq. (2) shows that the 
field at the detector drops off as the cube of the distance from the 
5 voxel of tissue. By applying a voxel magnetizing field that increases 
as the cube of the distance from the detector to the voxel, the field at 
the detector becomes independent of the distance from the voxet. In 
this case, the field is given by 

10 where a a and B Q are constants and ^ is the distance from the detector 
to the voxel. Thus, a magnetizing field that increases as the distance 
cubed from the detector "levels" the magnitude of signals from voxels 
at different locations relative to the detector. A three-dimensional 
magnetic susceptibility map can be generated from a three- 

15 dimensional map of the secondary magnetic flux of a magnetized 
phantom— such as a human (the corresponding space of detection is 
hereafter called the "sample space"). The magnetic susceptibility map 
is reconstituted using a Fourier transform algorithm. The algorithm is 
based on a closed-form solution of the inverse problem— solving the 

20 spatial distribution of magnetic dipoles from the extrinsic secondary 
flux from an array of magnetized voxels. The extrinsic magnetic flux 
is sampled at the Nyquist rate (twice the spatial frequency of the 
highest frequency of the Fourier transform of the magnetic 
susceptibility map of the phantom). 

25 In an embodiment of ihe present invention, nuclear magnetic 

resonance (NMR) is the means to measure the secondary magnetic 
field to provide the input 10 the magnetic susceptibility reconstruction 
algorithm. In this case, the measured secondary (RF) field is 
transverse to the magnetic flux including the local contribution due to 

30 the magnetic susceptibility of the voxel The RF field is detected in 
the near zone where r«X (or *r«|) t and the near fields arc quasi- 
stationary, oscillating harmonically as but otherwise static in 
character. Thus, the transverse RF magnetic field of each voxel is that 
of a dipolc, the maximum amplitude is given by Eq. (2) wherein the 

35 Larmor frequency of each voxel is shifted due to its magnetic 

susceptibility, and the magnetic moment along the z-axis. of Eq. 
(2) corresponds to the bulk magnetization M of each voxel. In terms 
of the coordinates of Eq. (2), an array of miniature RF antennas point 
samples the maximum dipolc component of the RF signal over the 

40 sample space such as the half space above (below) the object to be 
imaged wherein each RF signal is frequency-shifted by the 
perpendicularly oriented magnetic susceptibility moment of each 
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voxel. A three-dimensional magnetic susceptibility map is generated 
from a three-dimensional map of the secondary (RF) magnetic flux of 
a magnetized phantom— such as a human. The magnetic susceptibility 
of each voxel(s) is determined from the shift of the Larmor frequency 
5 due to the presence of the voxel(s) in the magnetizing field. The 
measurements of the spatial variations of the intensity of the 
transverse RF field of a given frequency is used to determine the 
coordinate location of each voxel(s). In one embodiment, the 
magnetic susceptibility map is reconstructed using a Fourier 

10 transform algorithm. The algorithm is based on a closed-form 

solution of the inverse problem— solving the spatial locations) of a 
magnetic dipole(s) of known magnetic moment (via the Larmor 
frequency) from the spatial variations of the extrinsic transverse 
secondary (RF) flux from the magnetized voxcl(s). The extrinsic (RF) 

13 magnetic flux is sampled at the Nyquist rate (twice the spatial 

frequency of the highest frequency of the Fourier transform of the 
magnetic susceptibility map of the phantom). In the limit, each 
volume element is reconstructed independently in parallel with all 
other volume elements such that the scan time is no sreater than the 
20 nuclear free induction decay (FID) time. 

The NMR type ReMSI scan performed on the object to be imaged 
including a human comprises the following steps: 

'The primary or magnetizing field is first determined over the voxels 
2d of the imaged volume (hereafter called "image space") in the 
absence of the object to be imaged. 

• The magnetic moments of nuclei including protons of the object to 

be imaged that arc aligned by the primary field are further aligned 
30 by a radio frequency (RF) pulse or series of pulses. . . ~ 

• The strength and duration of the rotating H, (RF) field that is 

resonant with the protons of the magnetized volume and is 
oriented perpendicularly io the direction of the magnetizing field is 
35 applied such that the final precession angle of the magnetization is 
90° (0 = 90*) such that the RF dipole is transverse to the primary 
magnetizing field and perpendicular to the RF maenetic field 
detector. 



40 



♦ The nuclei precess and emit RF energy at the frequency determined 
by the total magnetic flux at the voxel which is the superposition 
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of the primary (applied flux) and (he secondary (magnetic 
susceptibility flux). 

• The precessing nuclei undergo a free induction decay as they emit 
5 energy. 

• The magnetic intensity is time harmonic according to the magnetic 

flux at each voxel. The time dependent signals are Fourier 
transformed to give the spectrum. The magnetic moment 
10 corresponding to each frequency is determined (Eq. (21)). 

• The transverse RF Held is a dipole. The maximum intensity is 

determined for each frequency at each detector for at least one 
synchronized time point or short time interval to form a mat-ix for 

15 each frequency of the intensity as a function detector position in 
sample space. The measurements of the spatial variations of the 
transverse RF field of a given frequency is used to determine the 
coordinate location of each voxet(s). Thus, each matrix comprises 
the intensity variation over the sample space of the RF field of the 

20 bulk magnetization M of each voxel having the corresponding 
magnetic moment. To improve the signal to noise ratio of the 
image, this procedure may be repeated for multiple time points 
with reconstruction and averaging of the superimposed images, or 
the data may be averaged before reconstruction. 

25 

• The Fourier transform algorithm is performed on each frequency 
component over the detector array to map each bulk magnetization M 
corresponding to a magnetic moment to a spatial location or locations 
over the image space. The later case applies if the same magnetic 

30 moment is at more than one position in the image. The three- 
dimensional magnetic susceptibility map (the input function) is given 
as follows. With respect to the coordinate system of Eq. (2). (x, y t and 
z are the Cartesian coordinates. m_. the magnetic moment along the z- 
axis. of Eq. (2) corresponds to the bulk magnetization M of each voxel. 

35 and B is the magnetic flux due to the magnetic moment shown in 
FIGURE 9; the relationship to the NMR coordinate system is given in 
the Reconstruction Algorithm Section) the origin of the coordinate 
system, (0,0,0). is the center of the upper edge of the phantom. The 
phantom occupies the space below the plane x, y, z=Q (z<Q in the 

40 phantom space), and the sampling points lie above the plane (z>0 in 
the sampling space). The magnetic flux in the sampling space is given 
by multiplying the convolution of the input function with the system 
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function by the unitary function (one for ? >0 and zero for z <0). The 
input function can be solved in closed-form via the following 
operations: t? 

1- Measure the magnetic flux at discrete points in the sampling 
space. Each point is designated (x. y, z, flux) and each flux value is 
an element in matrix A. 

2. Discrete Fourier transform matrix A to obtain matrix B. 

3. Multiply each element (flux value) of matrix B by the 
corresponding inverse (reciprocal) value of the Fourier transform 
of the system function, Eq. (2), evaluated at the same frequency as 
the elemcm of the matrix A. This is matrix C. 

4. Generate matrix D by taking the discrete inverse Fourier 
transform of matrix C. 

5. Multiply each element of matrix D by the distance squared 
along the z-axis to which the element corresponds to generate the 
bulk magnetization M map. (This corrects the limitation of the 
sample space to z>0). 

6. In one embodiment, the voxels with a finite bulk magnetization 
25 M above a certain threshold or at an edge in the map are assigned 

the magnetic susceptibility (magnetic moment) determined from 
the shift of the Larmor frequency due to the presence of the 
voxel(s) in the uniform (nonuniform) magnetizing field. The other 
voxels arc assigned a zero value. This procedure is repeated for all 
30 Larmor frequencies. In the limit with sufficient Larmor frequency 
resolution, each volume element is reconstructed independently in 
parallel with all other volume elements such that the scan time is 
no greater than the nuclear free induction decay (FID) lime. 

35 7. In the case that the primary field is nonuniform, the magnetic 
. moment map is determined, then the magnetic susceptibility map 
of each magnetic moment is given by dividing the magnetic 
moment map by the primary magnetic field map on a voxel bv 
voxel basis and by subtracting one from each term (Eq. (23)). The 

40 total magnetic susceptibility map is the superposition of the 
separate maps for each magnetic moment corresponding to a 
unique Larmor frequency which is plotted and displayed. 
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BRIEF DESCRIPTION OF THE DRAWINGS 

FIGURE la is a schematic diagram of an apparatus according 
to the teach.ngs of the invention in which a magnet magnetizes a 
body to be tmaged, an RF generator excites the body, and an 
array of detectors detect radiation emitted by the body in 
response to the RF excitation radiation; 

FIGURE lb is a schematic diagram of an apparatus according to 
the teachings of the invention that employs a three-dimensional array 
of detectors to detect radiation emitted by nuclei of a body to be 
imaged m response to excitation of the nuclei by a radiation source- 
FIGURE I c is a plot of the Fourier transform H[k .k k ) of the * 
system function /,( w ) (Eq. (2)) corresponding to the z-component of 
a magnetic d.pole oriented in the e-dircction in accordance with the 
invention; 

FIGURE 2 is the plot of the field of a ring of dipoles of radius R 
and magnetic moment ,„ = 10' Gem" given by Eq. (1.14) as a function of 
radius R where the position of the center of the ring relative to the 
20 detector ,s the point (0,0.10) in accordance with the invention- 

FIGURE 3 is the plot of the field of a ring of dipoles of radius 
R - 0.2 cm and magnetic moment ,„ = )0 J Gem' given by Eq. (I 14) as a 
function of the distance between the detector at the origin and the 
center of the ring at the points (0,0. ; = 4 cm ,oz = 15 cm ) in accordance 
•i:> with the invention; 

FIGURE 4 is the plot of the field of a shell of dipoles of radius R 
and magneuc moment ,„ = | 0 ' Gcm> given by Eq. (11.17) as a function of 
radius R where the position of the center of the shell relative to the 
detector is the point (0.0. 10) in accordance with the invention. 

30 FIGURE 5 is the plot of the field of a shell of dipoles of radius 

R = 0.7cm and magnetic moment ,„ = |0 4 Con' given by Eq. (II 17) as a 
function of the distance between the detector at the oriain and the 
center of the shell at the points (0.0. : = 4 cm , ol = 15 cm) "in accordance 
with the invention; 

35 FIGURE 6 is the plot of the field of a sphere of dipoles of radius 

R and magnetic moment m = io J Gem' given by Eq (IV 16) as a 
function of radius R where the position of the center of the sphere 
relative to the detector is the point (0.0. 10) in accordance with the 
invention; 

40 FIGURE 7 is the plot of the field of a sphere of dipoles of radius 

R = 0.2 an and magnetic moment W = IO'G«n' given by Eq. (IV 16) as a 
function of the distance between the detector at the origin and the 
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center of ihe sphere at the points (0.0, ; = 4 cm,oz = 15cm) in 
accordance with the invention; 

FIGURE 8 shows a typical ReMSI (NMR) apparatus in accordance 
with the invention; 

5 FIGURE 9 shows the coordinate system (x. y. z) of Eq f>) with a 

primary field H„ and .he corresponding magnetic dipole both oriented 
parallel .o the z-axis wherein the z-component of the flux due to the 
z-onented d lp ole is measured at a detector according to Eq. (2) and 
shows the distances from the voxel ,o the detector in accordance with 
10 the invention: 

FIGURE 10 shows the stationary coordinate system of the ReMSI 
nuclear magnetic resonance (NMR) apparatus of FIGURE 8 
corresponds to the coordinate system of FIGURE 9 in accordance 
with the invention; 

15 FIGURE 11 shows the rotating NMR coordinate system (x R y 

z R ) and the stationary coordinate system (x, y. z) of the NMR detector 
corresponding to the coordinate system of FIGURE 9 and FIGURE 10 of 
a primary field H, oriented parallel to the z„-axis and the z-axis and 
the corresponding transverse RF magnetic d.pole oriented i„ the x K y - 

20 plane and periodically parallel to the y-axis wherein the spatial 

variation of the RF y-component of the flux due to the RF dipole is 
measured a. a detector according to Eq. (2) and shows the distances 

™ l ° ,hC dC,CCt0r in accorda "« with the invention; 

FIGURE 12 is a schematic of the three dimensional detector 
25 (antennae) array with respect to the stationary NMR coordinate 
system of FIGURE 10 which corresponds to the coordinate systems 
shown in FIGURES 9 and 11 in accordance with the invention, and 

FIGURE 13 shows the general process of reconstruction by 
reiteration; 

30 Further details regarding specific derivations, and calculations 

are provided ,n the attached appendices, wherein: 

APPENDIX I is the field produced by a rins of dipoles according 
to the present invention; " ° 

APPENDIX II is the derivation of the field produced by a shell of 
-55 dipoles according to the present invention; 

APPENDIX HI is the mathematical proof that the field produced 
by a shell of magnetic dipoles ,s different from that of a single dipole 
according to the present invention; 
f APPENDIX IV is the derivation of the field produced by a sphere 
40 of dipoles according to the present invention; 

APPENDIX V is ihe derivation of the Fourier transform of the 
system function of the z-componcnt of the magnetic field from a 
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dipole oriented in the direction of the z-axis used in the 
reconstruction process according to the present invention; 

APPENDIX VI is the derivation of S = HF®U(k z ) convolution from 
Eq. (44) used in a reconstruction process according to the present 
5 invention, and 

APPENDIX VII is the derivation of the Inverse Transform of Eq. 
(57) to Give Inverse Transform I. Eq. (58) T used in a reconstruction 
process according to the present invention. 

10 

DETAILED DESCRIPTION OF THE INVENTION 

An exemplary embodiment of an RcMSI apparatus 110 
according to the teachings of the present invention is shown in 

15 FIGURE la. The apparatus 110 includes a magnet 112. such as a 
superconducting magnet, that provides a primary or magnetizing 
field, to magnetize a body 114 to be imaged. A magnetometer 116 
determines the primary or magnetizing field in the volume to be 
occupied by the body t i.e., in the image space, in the absence of the 

20 body. One practice of the invention utilizes a magnetometer that 
employs NMR of protons in water for determining the primary 
radiation field at multiple points in the image space. A radiation 
source 118, such as a radio frequency generator, applies an RF pulse 
to the body in combination with RF antennae 118a when the body is 

25 placed in the magnetic field, to excite and thereby polarize selected 
nuclei of the body. The excited nuclei emit an RF radiation that a 
plurality of detectors 120, disposed in a plane above or below the 
object, can detect. The excitation pulse can be selected to rotate the 
magnetization of the nuclei, preferably by 90 degrees, with the 

30 respect to the primary field. In such a case, the RF radiation that the 
excited nuclei emit is primarily along a direction perpendicular to the 
plane of the detectors. The detectors can be selected to respond only 
to components of a magnetic field perpendicular to the plane in which 
they reside. Thus, such detectors can detect the emitted RF field 

35 without having interference from the primary magnetic field. A drive 
mechanism 122 moves the detectors 120 in a direction perpendicular 
to the plane of the detectors, to sample the external RF field over a 
three-dimensional volume. The separation of the detectors and the 
step size of the movement of the detectors along a direction 

40 perpendicular to the plane are selected such that the detectors sample 
the external RF field over a three-dimensional volume, i.e., the sample 
volume, at least at the Nyquist rate. Such a sampling advantageously 
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allows obtaining the spatial variations of the magnetic susceptibility 
of the body. 

A variety of detectors known in the art can be employed to 
5 construct the array of detectors 120. Such detectors include, but are 
not limited to, superconducting quantum interference devices 
(SQUID). Micromagnetic field sensors that employ galvanomagnetic 
effect due to the Lorentz force on charge carriers, semiconductor 
magnetic field sensors including MAGFET. magnetotransistors. Van 

10 dcr Pauw devices, integrated bulk Hall devices, such as the vertical 
MT and the lateral MT, silicon on sapphire and CMOS magnetodiodes, 
the magnetounijunction transistor, the carrier domain magnetometer, 
magnetic avalanche transistors, optoelectronic magnetic field sensors] 
and magnclorcsistivc magnetic field sensors. One preferred 

15 embodiment of the invention employs an array of SQUIDs, wherein 
each SQUID is cooled below the critical temperature after the 
magnetizing field is turned on and stabilized, to link the flux of the 
magnetizing field, and. thus, detect only the field generated by the 
body in response to an RF excitation. See Sarwinski, R. E., 

20 Superconducting Instruments. Cryogenics, Dec. 1977, pp.* 671-679, 
Another embodiment employs an impedance-matched array of RF 
antennas that are time-multiplexed to reduce cross talk among them. 

The process of magnetization of selected nuclei of the body can 

2 5 be better understood by referring to FIGURE 11 that shows a 

magnetizing field H 0 applied to a voxel 14a of the body 114. The RF 
excitation field is selected to be in a direction perpendicular to the 
field H fl and its magnitude is designated as 11,. NMR active nuclei of 
the voxel, such as protons, possess both angular momentum and a 

30 magnetic moment. Thus, the vector sum of the magnetic moments of 
all such NMR active nuclei present in the voxel I4a sive rise to the 
bulk magnetization of the voxel 14a. The bulk magnetization vector 
M (not shown) designates the collective contribution of a selected 
type of NMR active nuclei to the magnetization of the voxel 14a which 

35 corresponds to a RF parallel magnetic moment m. The RF excitation 
field is selected to be in resonance with the selected type of nuclei, 
e.g.. protons, of the voxel 1 4a, to rotate the bulk magnetization vector 
M. In a rotating frame representation, designated in FIGURE 11 as 
(x R , y R . z R ) and well known to those skilled^ in NMR r the magnetization 

40 vector M rotates about H, so long as the RF radiation is present. The 
rate of rotation of the magnetization vector about the applied field H, 
depends on the gyromagnetic ratio of the affected nuclei and the 
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magnitude of the excitation field. H,, The duration of die RF radiation 
can thus be selected to cause a rotation of the magnetization vector, 
initially aligned along H 0 , onto the x R y B -plane, for example by a 90° 
rotation. After the RF excitation field is turned off, the rotated 
5 magnetization M precesses about the local magnetic flux B (not 

shown) at the position of the voxel 14a, The local magnetic flux B is 
the sum of the primary magnetic flux ji 0 H. and the magnetic flux due 
to the magnetic susceptibility of the magnetized voxel, ^H, where // 0 
represents the permeability of free space, x is the magnetic 

10 susceptibility of the voxel, and H is the primary magnetic field 

strength at the position of the voxel The total magnetic flux is the 
sum of the contributions from the primary magnetic flux and the flux 
corresponding to the magnetic susceptibility of the magnetic voxel. 
The term (t+j) is herein defined as the total magnetic susceptibility. 

15 The magnetic susceptibility / causes a small perturbation of the 

magnetizing flux because the value of the magnetic susceptibility x »s 
typically 10" 4 to 10\ Thus, the magnetizing field musl be stable to at 
least to i part in I0 5 in order to produce a magnetic susceptibility 
image. Superconducting magnets arc particularly suited for providing 

20 such a stability because superconducting magnets typically provide a 
magnetic field stability of I part in 10* over a period of one month. 

The precession of the magnetization M of each voxel about the 
local flux B produces a radiating dipolar field corresponding to a 

25 magnetic moment m that has its maximum intensity along the y R 

direction. The precession frequency of the magnetic moment called 
the Larmor frequency is dependent on the local magnetic flux which 
has a contribution corresponding to the magnetic susceptibility of the 
voxel (See Eq. (21) and Eq. (22)). A detector 20a can.be, for example, 

30 an RF antenna pointing along the y direction to respond selectively to 
the radiating RF field from all of the voxels. A Fourier transform of 
this external RF field recorded as a function of lime produces its 
frequency components ("Larmor frequencies"). Each frequency 
component can provide the total local magnetic flux in voxels 

35 producing that frequency component. The total local magnetic flux in 
turn is related to the susceptibility of the material within these voxels 
(See Eq. (22)). Thus, the values of magnetic susceptibility can be 
gleamed from frequency components of the external RF radiation. 
The next step in the process of determining the magnetic 

40 susceptibility image involves determining locations of the voxels 
contributing to each frequency component. The variation of the 
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maximum intensity of the RF field at a particular frequency over the 
three-dimensional sample space is used to determine the locations of 
the voxels contributing to radiation at that frequency. Preferably, the 
magnetic susceptibility of each voxel is reconstructed independently 
5 in parallel with all other voxels such that the scan time is no greater 
the time for the excited nuclei to return to their unexcitcd state called 
the nuclear free induction decay (FID) time. This may be achieved by 
using sufficiently high resolution of the Larmor frequencies such that 
each frequency component corresponds to a sinele voxel 

10 

With further reference to FIGURE la, an analogue to digital 
converter CA/D converter") 124 converts the analogue signal outputs 
of the detectors 120 into digital signals. A preferred embodiment of 
the invention utilizes at least a 12 bit A/D converter to digitize the 

15 output signals. A reconstruction processor 126 receives the digital 
signals, and determines the magnetic susceptibility map of the body. 
The reconstruction processor 126 includes a Fourier transform 
processor 126a that obtains the frequency components of the external 
RF field. In addition, a spatial locator 126b. which is a part of the 

20 reconstruction processor 126. employs the variation of the maximum 
intensity of the external RF field over the sample space (ihc three 
dimensional space sampled by the detectors) to locate the positions of 
the voxels in the image space producing a particular frequency 
component of the external RF field, in a manner described in detail 

25 below (See the Reconstruction Algorithm Section). An display 

processor 128 displays a two-dimensional or a three-dimensional 
image corresponding to the variations of the magnetic susceptibility 
of the body. 

30 Thc spatial locator 126b employs an algorithm, described in 

detail below in the Reconstruction Algorithm Section, to determine the 
positions of the voxels of the magnetized body that produce each of 
the frequency components of the external RF field. It is demonstrated 
in the Uniqueness of Reconstruction Section that a magnetic field 

35 produced by a geometric distribution of dipoles is unique. Therefore, 
a unique spatial distribution of magnetic dipoles. such as those 
corresponding to the bulk magnetization of each voxel due to 
precessing nuclei, gives rise to a unique magnetic field. Thus, the 
measured external RF field can provide a unique solution for the 

40 spatial distribution of magnetic dipoles. i.e.. magnetized voxels 
comprising excited nuclei, in the body. 
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Referring again to FIGURE la, the spatial locator 126b preferably employs 
a Fourier Transform Reconstruction Algorithm to map a bulk magnetic moment m 
corresponding to a particular frequency component of the external RF field onto a 
spatial location or locations in the body. The latter case applies if more than one 
5 location in the body gives rise to a particular frequency component of the external 
RF radiation. 

According to a preferred embodiment of the invention, the system 110 
employs an open design magnet, such as a Helmholtz coil, to allow positioning of 
10 the detector array 120 sufficiently close to the body 1 14 to be imaged. The 
orientation of the magnetic field with respect to the body can be selected to 
optimize the signal to noise ratio of the signals detected by the array of detectors 
1 20. For example, in the case of imaging a patient body, the primary magnetic 
field can be selected to be coaxial with the body, or it can alternatively be 

1 5 perpendicular to the body axis. 

The ReMSI apparatus of FIGURE la provides a number of advantages. In 
particular, because the magnetic susceptibility of tissue is sensitive to the 
concentration of oxyhemoglobin relative to deoxyhemoglobin in tissue, the 
apparatus of the invention is particularly suited for imaging cardiopulmonary and 

20 vascular systems. Further, the apparatus of the invention can achieve a typical 
resolution of 10 * cm 3 , thereby permitting physicians to view human anatomy and 
pathology in a manner not available with conventional imaging techniques. In 
addition, the fast data acquisition possible with the ReMSI technique reduces 
motion artifacts in the image. Another important advantage of the present 

25 technique is that it is non-invasive. In contrast, many of the conventional 

techniques, such as X-ray angiography, require surgery. Many conventional non- 
invasive imaging techniques, such as nuclear medicine and magnetic resonance 
angiography, provide at least an order of magnitude lower resolution than the 
technique of the invention, and are much more susceptible to motion artifacts, 

30 further, the present technique provides a three dimensional image that can be 
displayed from any perspective. 

An alternative embodiment of the ReMSI apparatus of the invention, 
shown in FIGURE lb, employs a three-dimensional array of detectors 230, 
spaced apart to detect spatial variations of the emitted RF radiation at least at 

35 the Nyquist frequency. A magnet 2 1 2 provides a magnetizing field in a 
volume to be occupied by a body, i.e., image space. A magnetometer 216 
measures the magnetizing field at a plurality of positions in the image space 
in the absence of the body. As in embodiment of FIGURE la, the body 214 
to be imaged is placed in the magnetizing field provided by the magnet 212. 

40 An RF generator 2 1 8 in combination with an RF antenna 2 1 8a apply an RF 
pulse to the body to polarize selected nuclei of the body. The polarized 
nuclei emit an RF radiation. The three-dimensional array of detectors 230 
provide output signals in response to the RF radiation emitted 
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RF radiation al least at the Nyquist frequency. A magnet 212 provides a 
magnetizing field in a volume to be occupied by a body ? i.e., image space. A 
magnetometer 216 measures the magnetizing field at a plurality of positions in the 
image space in the absence of ihc body. As in embodiment of FIGURE la, the 
5 body 2 1 4 to be imaged is placed in the magnetizing field provided by the magnet 
212. An RF generator 2 1 8 in combination with an RF antenna 2 1 8a apply an RF 
pulse to the body to polarize selected nuclei of the body. The polarized nuclei emit 
an RF radiation. The three-dimensional array of detectors 230 provide output 
signals in response to the RF radiation emitted by the body. A digitizer 224 

1 0 digitizes the output signals and sends the digital signals to a construction processor 
226 that determines variations of the magnetic susceptibility of the body in a 
manner similar to that described in connection with the embodiment of FIGURE la. 
A display processor 228 receives the information regarding the spatial variations of 
the magnetic susceptibility from the construction processor 226. and provides a 

15 two-dimensional or a three-dimensional magnetic susceptibility image of the body. 
Employing a three-dimensional array of detectors advantageously decreases the 
acquisition time because the emitted RF signal over the entire three-dimensional 
sample space is detected at once. The shortening of the acquisition lime in turn 
reduces motion artifacts in the susceptibility image. 

20 Uniqueness of Reconstruction 

The reconstruction processor 126 of FIGURE 1 constructs the magnetic 
susceptibility map according to a defined sequence and program. Those of 
ordinary skill will be able to properly program the processor in light of the 
following description. Those of ordinary skill will also recognize that the steps to 

25 reconstruct the magnetic susceptibility map can be implemented in either hardware 
or software. All matter has a permeability different from that of free space that 
distorts an applied magnetic flux. This property is called magnetic susceptibility. 
An object, herein called a phantom, can be considered as a collection of small 
volume elements or voxels. When a magnetic flux is applied to the phantom, each 

30 voxel generates a secondary magnetic flux external to the phantom. The strength 
of the secondary magnetic flux varies according to the strength of the applied flux, 
the magnetic susceptibility of the material within the voxel and the distance of the 
external location relative to the voxel. The net magnetic flux at a point extrinsic to 
the phantom is a sum of the applied flux and the contributions of each of the voxels 

35 within the object. This flux is point sampled over a three dimensional space, and 
the magnetic susceptibility of each voxel is solved by a reconstruction algorithm. 
The resulting image is called the magnetic susceptibility map. 

The secondary magnetic field due to magnetized tissue has to be 
40 modeled as noninteracting dipoles aligned with the imposed field. It 
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is demonstrated below that the field of any geometric distribution of 
dipoles is unique, and the superposition principle holds for magnetic 
fields; therefore, a unique spatial distribution of dipoies gives rise to a 
unique secondary magnetic field, and it is further demonstrated 
5 below that this secondary field can be used to solve for the magnetic 
susceptibility map in closed form. It follows thai this map is a unique 
solution. To prove that any geometric distribution of dipoles has a 
unique field, it must be demonstrated that the field produced by a 
dipole can serve as a mathematical basis for any distribution of 
10 dipoles. This is equivalent to proving that no geometric distribution 
of dipoies can produce a field which is identical to the field of a 
dipole. 

By symmetry considerations, only three distributions of uniform 
dipoles need to be considered. A magnetic dipole has a field that is 

15 cylindrically symmetrical. A ring, a shell, a cylinder, and a sphere of 
dipoles are the only cases which have this symmetry. A cylinder is a 
linear combination of rings. Thus, the uniqueness of the dipole field is 
demonstrated by showing that it is different from that of a ring, a 
shell, and a sphere. The uniqueness of the dipole from the cases of a 

20 ring, a shell, and a sphere of dipoles is demonstrated in APPENDIX I, 
APPENDIX IL and Appendix IV, respectively. The plot of the three 
cases of the field of a ring> shell, and a sphere of dipoles each of 
radius R and magnetic moment m=\0 4 Gem* given by Eq. (1.14). Eq. 
(H-17), and Eq. (IV.16) of APPENDIX I t APPENDIX II. and APPENDIX 

25 IV as a function of radius R where the position of the center of each 
distribution relative to the detector is the point (0.0,10) is given in 
FIGURES 2. 4. and 6. respeciively. Since the fields vary as a function 
of radius R . the dipole field is not equivalent to these distributions of 
dipoles. k is further mathematically proven in APPENDIX III that the 

30 field produced by a shell of magnetic dipoles is different from that of 
a single dipole. AH other fields are a linear combination of dipoies. 
Thus, the dipole is a basis element for the reconstruction of a 
magnetic susceptibility image. Since each dipole to be mapped gives 
rise to a unique field and since the total field at a detector is the 

35 superposition of the individual unique dipole fields, linear 

independence is assured; therefore, the ReMSI map or image is 
unique. In other words, there is only one solution of the ReMSI image 
for a given set of detector values which spatially measure the 
superposition of the unique fields of the dipoles. This map can be 

40 reconstructed using the algorithms described in the Reconstruction 
Algorithm Section. 
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Reconsiruciion of the Magnetic Susceptibility Map 

The volume to be imaged is divided into volume elements called 
voxels, and the magnetized voxel 302 shown in FIGURE 9 with 
magnetic moment density X B is modeled as a magnetic dipole where 
X is the average volume magnetic susceptibility of the voxel and B is 
the magnetizing flux at the position of the dipole corresponding to a 
magnetizing field H„. Data comprising the z-componcnt of the 
magnetic field of a dipole oriented in the z-dircction is acquired by 
detectors 301 in the three dimensional sample space comprising the 
xy-plane and the positive z-axis as shown in FIGURE 9. The magnetic 
moment. w : . of each voxel whhin the phantom is given by the product 
of its volume magnetic susceptibility, x , the magnetizing flux oriented 
along the z-axis. B. r and the volume of the voxel. AV. 

The magnetic moment m of each voxel is a magnetic dipole. And. the 
phantom can be considered to be a three-dimensional array of 
magnetic dipoles. At any point extrinsic to the phantom, the z- 
component of the secondary flux. B, from any single voxel is 

20 ti =m _J^ziLz^L fls , 

where x. y, and z are the distances from the center of the voxel to the 
sampling point. It is shown in APPENDICES I-IV that no geometric 
distribution of magnetic dipoles can give rise to Eq. (2). Therefore, 
the flux of each magnetic dipole (voxel contribution) forms a basis set 

2o for the flux of the array of dipoles which comprise the magnetic 
susceptibility map of the phantom. 

Eq. (2) is a system function which gives the magnetic flux 
output in response to a magnetic dipole input at the origin. The 
phantom is an array of spatially advanced and delayed dipoles 

30 weighted according to the magnct.c susceptibility of each voxel; this is 
the input function. The secondary flux is the superposition of 
spatially advanced and delayed flux, according to Eq. (2); this is the 
output function. Thus, the response of space to a magnetized 
phantom is given by the convolution of Eq. (2) with the series of 

35 weighted, spatially advanced and delayed dipoles representing the 
susceptibility map of the phantom. The discrete signals arc recorded 
by a detector array over the sample space comprising the xy-plane 
and the positive z-axis. ' * 

In an embodiment of the present invention, nuclear magnetic 

40 resonance (NMR) is the means to measure the secondary maanetic 

field to provide the input to the magnetic susceptibility reconstruction 
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algorithm. In this case, the measured secondary (RF) field is 
transverse to the magnetic flux including the local contnbution due to 
the magnetic susceptibility of the voxel. The RF field is detected in 
the near zone where r<<A (or kr«\), and the near fields are quasi- 
~* stat,ona 'y. oscillating harmonically as e' iu , but otherwise static in 

character. Thus, the transverse RF magnetic field of each voxel is that 
of a dtpole, the maximum amplitude- is given by Eq. (2) wherein the 
Larmor frequency of each voxel is shifted due to its magnetic 
susceptibility, and the magnetic moment along the z-axis of Eq 

10 (2) corresponds to the bulk magnetization M of each voxel The 

reconstruction may be multiplexed by using the NMR signals as the 
input The NMR frequency is determined by the magnetizinc flux 
which includes a contribution from the magnetic susceptibility of the 
voxel Each voxel or voxels which have a unique local magnetic flux 

15 w,H have a unique a NMR frequency. The spatial varia.ion of the 
intensity of the transverse RF signal of this unique frequency is the 
input for the reconstruction algorithm. This procedure is performed 
either in parallel or series for each data set of signals of each unique 
frequency. The magnetic susceptibiluy image is the superposition of 
20 the independent images, each of which corresponds to a unique 
frequency. 

Frequency Multiplexing By Nuclear Magnetic Resonance Measurement 
25 of thc Secondary Field Over Hie Sample Space 

The NMR active nuclei including protons posses both angular 
momentum and a magnetic moment. When nuclei are placed in a 
static magnetic field H 0 . they precess about ihc field at a frequency 
proportional to the magnitude of H 0 . The bulk magnetization M of 

30 each voxel comprises the vector sum of the magnetic moments from 
all of the nuclei in each voxel. If the prccessing nuclei arc then 
subjected to an additional rotating (RF) field H„ which is synchronous 
with the precession, their magnetic moments and thus M will precess 
about H, and rotate away from the primary field H, bv an angle * in 

35 a coordinate frame which rotates at the Larmor frequency (31 The 
precession about H, continues as long as H, exists. The final value of 
P then depends on the strength of H,. which determines Ihc 
precession rate, and the time for which it is turned on. The nuclei 
absorb energy as they change their orientation. This is known as 
40 nuc car magnetic resonance (NMR). The temperature of the nuclei or 
nuclear spin system rises during absorption of energy. When the H 
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Held is removed, the spin system cools down until it is thermal 
equilibrium with its environment. The exponential relaxation of the 
spin system temperature to that of the surrounding lattice is called 
spin-lattice relaxation and has a time constant T, where a time 
constant is deftned as the time it takes for 63% of the relaxation to 
occur. The NMR signal may also decay because the nuclei initially in 
phase following the H, pulse get out of alignment with each other or 
dephase by local interactions with the maenetic fields of neighbor 
nuclei. The dephasing of the NMR signal is due to differing precession 
rates effected by the Jocal interactions and is described by an 
exponential time constant T r The main source of NMR image ( also 
called magnetic resonance images (MR[)) contrast is T, and T. which 
depend on tissue types. In an embodiment of the present invention. 
T, and T : are used to provide image contrast as described in the 
15 Contrast Section. 

In other words, when nuclei are placed in a static magnetic field 
H 0 and then subjected to an additional rotating (RF) field H,, which is 
synchronous with their precession. M will preccss about H, and 
rotate away from the primary field H 0 by an angle The magnitude 
20 of M is a maximum initially and decays with time. This occurs by 

emission of the same multipolarity radiation that it absorbed and by 
transfer of energy to the surrounding lattice. The intensiiv of the 
radiation is a function of M and the coordinate position re'lative to the 
RF emitting voxel. In the present invention, (he measurement of the 
2^ intensity of the RF signal is performed over time and space. The 
signal as a function of lime at a given detector position is Fourier 
transformed to give the spectrum. The magnetic susceptibilitv of 
each voxel is determined from the shift of the Lanr.or frequency due 
to the presence of the voxel in the magnetizing field. The intensity 
30 vanation of the transverse RF field over space is used to determine 
the coordinate location of each voxel. 

The nature of the RF field can be determined from Maxwell's 
equations applied to a sinusoidal current. With a sinusoidal current 
J(x') confined to small region compared with a wavelength, the 
35 solution of the vector potential A(x) is [4] 



A(x)= ^ J(x,) irrxi rfV <'<>> 

where k = j is the wavenumber. and a sinusoidal time dependence is 

understood. The magnetic induction is siven by 

B = VxA (17) 
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while, outside the source, the electric field is 

E 4 VXB (18) 
For a source of dimension d + the fields in the far zone defined by 
d«X«r are transverse to the radius vector and fall off as r~\ 
5 typical of radiation fields. For the near zone where r«X (or *r«l), 
the exponential in Eq. (16) can be replaced by unity. Then the vector 
potential is given by 

lim x I v 4 it Y. (6,6) { . . 

l^To A{x) - cZtt-^IW'Ki'.W* ( 1 9) 

This shows that the near fields are quasi-stationary, oscillating 

10 harmonically as e~ ,a \ but otherwise static in character. 

Nuclear magnetic resonance (NMR), which is commonly called 
magnetic resonance imaging (MRI), is a means to measure the 
primary and secondary magnetic fields to provide the input to the 
magnetic susceptibility reconstruction algorithm. The proton 

15 gyromagnetic ratio y F /2n is 

y p 1 2/r = 42.57602 MHz T % (20) 
The NMR frequency / is the product of the proton gyromagnetic ratio 
given by Eq. (20) and the magnetic flux B. 

/ = y, / 2/rB = 42.57602 MHz T~ l B (21) 

20 A typical flux for a superconducting NMR imaging magnet is 0.25 T. 
According to Eq. (21) this corresponds to a radio frequency (RF) of 
10.6 MHz which corresponds to a wavelength of 28.3 m. In the present 
invention, each RF antennae of an array is located at a distance of 
about iOcm from the voxels within the image space. Thus, the RF 

25 field is detected in the near zone where r«X (or *r«l), and the 
near fields according to Eq. (19) arc quasi-stationary, oscillating 
harmonically as but otherwise static in character. The transverse 

RF magnetic field of each voxel is that of a RF dipole, the maximum 
amplitude is given by Eq. (2) wherein the Larmor frequency of each 

30 voxel is shifted due to its magnetic susceptibility, and #w. % the 

magnetic moment along the z-axis, of Eq. (2) corresponds to the bulk 
magnetization M of each voxel. 

The magnetic susceptibility gives rise to a magnetic moment 
and a secondary magnetic field at the position of the voxel as well as 

35 external to the object to be imaged. In fact, magnetic susceptibility 

has been imaged by conventional MRI under extraordinary conditions 
[5j. In the case of NMR detection of this secondary field according to 
the present invention, the transverse RF field is measured. An array 
of miniature RF antennas point samples the dipole component (Eq. 
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(2)) of the RF signal due to the perpendicularly oriented maenctic 
susceptibility moment of each voxel over the sample space such as 
the half space above the object to be imaged. The magnetic 
susceptibility of each voxel is determined from the shift of the 
5 Larmor frequency due to the presence of the voxel in the magnetizing 
field. The intensity variation of the transverse RF field over space is* 
used to determine the coordinate location of each voxel. The RF field 
is the near field which is a dipole (with a maximum amplitude given 
by Eq. (2)) that serves as a basis clement to form a unique 
10 reconstruction. The geometric system function corresponding to a 
dipole which determines the spatial intensity variations of the RF 
field is a band-pass for * p = A ; . 1„ the limit, each volume element is 
reconstructed independently in parallel with all other volume 
elements such that the scan time is no greater than the nuclear free 

15 induction decay (FID) time. 

In the case of NMR detection, the RF frequency signal from a 
magnetized object (Eq. (21)) is due to the total (lux B which is the 
sum of the flux from the primary magnetizing field /r 0 H and the 
secondary flux t i oX ll from the magnetic susceptibility x of the voxel- 

20 B = p 0 (| +/r ) H {22) 

where //„ is the permeability of free space and H is the primary 
magnetizing field strength. Eq. (22) expressed in the form of Eq. (1) is 

m. =(I + *)BAV ( 2 3) 

where ,i 0 H is the magnetizing flux and is herein defined as the 

25 total magnetic susceptibility. It can be seen from Eq. (23) that the 
total magnetic susceptibility map can be determined from the total 
magnetic moment map if the magnetizing flux map is known. Since 
* = i0~* -10 * «l. the magnetic susceptibility -/ causes a small 
perturbation to the magnetizing flux. Thus, the magnetizing flux must 

30 be stable to at least the magn.tude of the signal contribution due to 
the magnetic susceptibility. A feature of superconducting magnets is 
their extreme stability. The object to be imaged is magnetized' with a 
highly stable magnet such as a superconducting magnet." In this case, 
a magnetic field stability of lO'* % over a month's time is feasible. It 

35 can further be seen from Eq. (21) and Eq. (23) that the total magnetic 
moment map can^be determined from the NMR frequencies as input. 
Since % = I0 _J -I0- S «l ! the magnetic susceptibility % causes a small 
perturbation to each NMR frequency (Eqs. (21) and (22)). However, 
extremely high frequency resolution is routine. The frequency 

40 resolution of RF signals exceeds 10"* %. Small antennas measure the 
RF signals as point samples without significant decrease in the signal 
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to noise ratio relative to large antennas by using impedance matching 
while minimizing resistive losses by using superconducting reactance 
elements, for example. In an embodiment, cross talk between 
antennas is ameliorated or eliminated by time multiplexing the signal 
5 detection over the array of antennas. 

In an embodiment, the strength and duration of the rotating H ( 
(RF) field that is resonant with the protons of the magnetized volume 
and is oriented perpendicularly to the direction of the magnetizing 
field is applied such that the final precession angle of the 

10 magnetization is 90° (# = 90°) such that the RF dipole is transverse to 
the static magnetizing field and perpendicular to the RF magnetic field 
detector 20a as shown in FIGURE 11, The transverse RF field is a 
dipole. The intensity is determined for each frequency at each 
detector for at least one synchronized (common) time period or short 

15 time interval (the period corresponding to the Larmor frequency is 
KT'-KT'sec compared to the 7, (spin/lattice) and 7\ (spin/spin) 
relaxation times of 1-5 sec) to form a matrix for each frequency of the 
intensity as a function detector position in sample space (i.e. spatial 
variation of the intensity at a determined common lime point). The 

20 magnetic moment corresponding to each frequency is given by Eq. 
(21). Thus, each matrix comprises the intensity variation over the 
sample space of the RF field of the bulk magnetization M of each 
voxel having the corresponding magnetic moment. The rotating NMR 
coordinates and the stationary NMR coordinates are both shown in 

25 FIGURE 1L According to Eq. (2), the signal as a function of time which 
is Fourier transformed arises from each transverse RF magnetic dipole 
oriented in the x R y R -plane (rotating NMR coordinates) which is 
periodically parallel to the y-axis and rotates about the primary field 
■I, oriented parallel to the z-axis at the Larmor frequency (stationary 

30 NMR coordinates). For the stationary NMR case, the y-axis 

corresponds to the z-axis of Eq. (2). and the magnetic moment 
along the z-axis, of Eq. (2) corresponds to the bulk magnetization M of 
each voxel. 

In addition, the spectrum comprises an additional component 
35 corresponding to the measurement of the RF magnetic dipole which is 
periodically parallel to the x-axis and rotates about the primary field 
H 9 oriented parallel to the z-axis at the Larmor frequency wherein 
the x-axis of the stationary NMR coordinates corresponds to the x-axis 
of Eq. (2). The components rotate at the same frequency but are 
40 orthogonal. The orthogonality corresponds to a phase angle of the 
Fourier transform of 90°; thus, the spectrum contains a real and an 
imaginary part at the each Larmor frequency. The magnitude of the 
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amplitude of the signal due to the dipole oriented parallel to the v 
axis (z-axis of Eq. (2)) exceeds that of the dipole oriented parallel to 
the x-ax.s (x-axis of Eq. (2)). In a „ embodiment, the spatial variation 
over the detector array of the component at each Larmor frequency 
with the maximum amplitude is used to determine the coordinate 
locauon of each voxel using the geometric system function of the 
detection of the 2-componeni of a z-oriemed dipole as given in the 
Reconstruction Algorithm Section. This geometric system function is . 
band-pass for *, = *.. I„ another embodiment, the spatial v ar ,ation 
over « he detector array of the lesser magnitude orthogonal component 
at each Larmor frequency , s used to determine the coordinate locauon 
of each voxel using the corresponding geometric system function as 
described in the Reconstruction Algorithm Section. The system 
funct.on for the case of the minor orthogonal component in terms of 
the coordinates given ,n FIGURE 9. corresponds to the detection of the 
x-component of a z-oriented dipole. The geometric system function is 
given by the function of the x-component of the flux B from a z- 

sm g r d voxe P i° ,C " P ° im "'""^ '° ,,>C imagC SpaCC - fr ° m an y 



B3.rz 
= m 



U' + v'+r) 1 '- (24) 
where x, y, and z are the distances from the center of the voxel to the 
sampling point and m . is the magnetic moment along the z-axis 
corresponding to the bulk magnetization M of each voxel for the NMR 
esse . 

The magnetic susceptibility map is obtained by first 
reconstructing ihe magnetic moment map. The ReMSI scan 

fhTfolting sly*" 1 10 bC ima§ed 3 hUma " 

♦ The primary or magnetizing field is first determined over the imaec 
space in the absence of the object to be imaged In one 
embodiment. NMR of a proton containing homogeneous phantom of 
known susceptibility such as water is determined on a point by 
point bas.s to map the primary field, or the magnetic flux at 
multiple points is obtained simultaneously using the reconstruction 
algorithm described herein via determining the magnetic moment 
map over all or a portion of the image space, and then determinine 
the field us.n- the known magnetic susceptibility of the phantom 
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• The magnetic moments of nuclei including protons of the object to 

be imaged that are aligned by the primary field arc further aligned 
by a radio frequency (RF) pulse or series of pulses. 

5 ■ The strength and duration of the rotating II, (RF) field that is 
resonant with the protons of the magnetized volume and is 
oriented perpendicularly to the direction of the magnetizing field is 
applied such that the final precession angle of the magnetization is 
90 6 (<p = 90°) such that the RF dipole is transverse to the primary 
10 magnetizing field and perpendicular to the RF magnetic field 
detector. 

• The nuclei precess and emit RF energy at the frequency determined 

by the total magnetic flux at the voxel which is the superposition 
15 of the primary (applied flux) and the secondary (magnetic 
susceptibility flux). 

• The precessing nuclei undergo a free induction decay as they emit 

energy. 

20 

• The magnetic intensity is time harmonic according to the magnetic 

flux at each voxel. The time dependent signals are Fourier 
transformed to give the spectrum. The magnetic moment 
corresponding to each frequency is determined (Eq. (21)). 

25 

• The transverse RF field is a dipole. The maximum intensity is 

determined for each frequency at each detector for at least one 
synchronized (common) time period or short time interval to form 
a matrix for each frequency of the intensity as a function detector 

30 position in sample space (i.e. spatial variation of the intensity at a 
determined common time point). Thus, each matrix comprises the 
intensity variation over the sample space of the RF field of the 
bulk magnetization M of each voxel having the corresponding 
magnetic moment. (Typically, if more than one voxel has the same 

35 magnetic susceptibility and the magnetization field is uniform, 

then they will have essentially an equivalent bulk magnetization 
M. In this case and the case that only one voxel gives rise to a 
given Larmor frequency, the matrix comprising the intensity 
variation over the sample space of the RF field of the bulk 

40 magnetization M of each voxel having the corresponding magnetic 
moment is essentially the matrix comprising the intensity variation 
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over the sample space of the secondary field due to the magnetic 
susceptibility with each term multiplied by a constant.) 

* The Fourier transform algorithm given in the Fourier Transform 

5 Reconstruction Algorithm Section is performed on each frequency 
component over the detector array to map the bulk magnetization 
M of each voxel having the corresponding magnetic moment to a 
spatial location or locations over the image space. The later case 
applies if the same magnetic moment is at more than one position 
10 in the image. 

• In one embodiment, the voxels with a finite bulk magnetization M 

above a certain threshold or at an edge in the map are assigned the 
magnetic susceptibility (magnetic moment) determined from the 

15 shift of the Larrnor frequency due to the presence of the voxel(s) 
in the uniform (nonuniform) magnetizing field. The other voxels 
are assigned a zero value. This procedure is repeated for all 
Larrnor frequencies. In the limit with sufficient Larrnor frequency 
resolution, each volume element is reconstructed independently in 

20 parallel with all other volume elements such that the scan lime is 
no greater than the nuclear free induction decay (FID) time. 

♦ In the case that the primary field is nonuniform, the magnetic 

moment map is determined, then the magnetic susceptibility map 
25 of each magnetic moment is given by dividing the magnetic 

moment map by the primary magnetic field map on a voxel by 
voxel basis and by subtracting one from each term (Eq. (23)). 

♦ The total magnetic susceptibility map is the superposition of the 
.30 separate maps for each magnetic moment corresponding to a 

unique Larrnor frequency which is plotted and displayed. 

* The quality of the image (i.e. the signal to noise ratio of the image) 

can be increased by repeating the reconstruction over multiple 
35 time points wherein each data set of a given time point is the set of 
matrices of the intensity variation over ihe sample space of the RF 
field of the bulk magnetization M of each voxel having the 
corresponding magnetic moment at that synchronized (common) 
time point. 

40 

System Construction 
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The ReMSI system 110 of FIGURE la and the system of FIGURE 8 includes (1) 
a magnet including a superconducting magnet to magnetize a volume of an object 
or tissue to be imaged, (2) a magnetometer to determine the primary or 
magnetizing flux over the image space in the absence of the object to be imaged, 
5 (3) a two or three dimensional detector array which is translated over the third 
dimension or a three dimensional detector array wherein each detector is 
selectively responsive to a single component of the secondary magnetic field 
which is produced by the magnetized voxels of the object including tissue where 
the component of the field to which the detector is responsive determines the 
10 geometric system function which is used in the reconstruction algorithm 
discussed in the Reconstruction Algorithm 

Section, and the detectors ideally point sample the secondary magnetic field at the 
Nyquist rate over the spatial dimensions which uniquely determine the 
susceptibility map which is reconstructed from the measurements, 4.) an analogue 

1 5 to digital converter to digitize the secondary magnetic field measurements, 5.) a 
reconstruction algorithm processor including a Fourier transform processor to 
map the magnetic moments to the image coordinates and convert each magnetic 
moment map into a magnetic susceptibility map, and 6.) image processors and a 
display such that the ReMSI image can be rotated in space to be displayed from 

20 any perspective as a three dimensional or two dimensional (tomographic) image. 

Magnetizing Field 
In an embodiment of the present invention, the applied magnetizing field 
which permeates the object to be imaged including tissue is confined only to that 
region which is to be imaged. The confined field limits the source of signal only 

25 to the volume of interest; thus, the volume to be reconstructed is limited to the 
magnetized volume which sets a limit to the computation required, and eliminates 
end effects of signal originating outside of the edges of the detector array. In the 
NMR case, the field having a steep gradient at the edges limits the imaged region 
by providing a range of Larmor frequencies wherein the data is rejected at a 

30 threshold Larmor frequency. 

A magnetic field gradient can be applied in the direction perpendicular to 
the plane of the detector array described below to alter the dynamic range of the 
detected signal, as described in the Altering the Dynamic Range Provided by the 
System Function Section or in other orientations to improve image quality. The 
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magnetizing means can also possess a means to add a component of 
modulation to the magnetizing field at frequencies below those which 
would induce eddy currents in the tissue. Such modulation would 
cause an in-phase modulation of the RF secondary magnetic field 
5 signal which can be filtered to improve the signal to noise ratio, or the 
modulation over space can be used to encode spatial information. 

Dynamic Range 

10 As described in the Nyquist Theorem with the Determination of 

the Spatial Resolution Section, the system function which effects a 
dependence of the signal on the inverse of the separation distance 
between a detector and a dipole cubed band-passes the magnetic 
susceptibility function; however, the resulting signals are of large 

15 dynamic range requiring at least a 12-bit AJD converter. A 12 bit 
A/D converter is sufficient with a magnetizing design which exploits 
the dependence of the signal on the strength of the local magnetizing 
field. In an embodiment of the case shown in FIGURE 9, a quadratic 
magnetizing field along the z-axis is applied in the z-dircction such 

20 that the fall-off of the signal of a dipole with distance is compensated 
by an increase in the signal strength due to a quadratically increasing 
local magnetizing field. The reconstruction algorithm is as discussed 
in the Reconstruction Algorithm Section; however, each element of the 
matrix containing the magnetic moment values is divided by the 

25 corresponding clement of the matrix containing the values of the flux 
magnetizing the corresponding magnetic susceptibility. These 
operations are described in detail in the Altering the Dynamic Range 
Provided by the System Function Section. Eq. (2) shows that the field 
at the detector drops off as the cube of the distance from the voxel of 

30 tissue. By applying a voxel magnetizing field that increases as the 
cube of the distance from the detector to the voxel, the field at the 
detector becomes independent of the distance from the voxel. In this 
case, the field is given by 

B ;=^>(«- + ^P (25) 

35 where a m and B 0 arc constants and z n is the distance from the detector 
to the voxel. Thus, a magnetizing field that increases as the distance 
cubed from the detector "levels" the magnitude of signals from voxels 
at different locations relative to the detector. 

In the NMR case shown in FIGURE 11, a quadratic magnetizing 

40 field along the y-axis is applied in the z-dircction such that ~thc fall- 
off of the signal of a dipole with distance is compensated by an 
increase in the signal strength due to a quadratically increasing local 
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magnetizing field. The reconstruction algorithm is as discussed in the 
Reconstruction Algorithm Section; however, each element of the 
matrix containing the magnetic moment values is divided by the 
corresponding element of the matrix containing the values of the flux 
5 magnetizing the corresponding magnetic susceptibility. These 

operations arc described in detail in the Altering the Dynamic Range 
Provided by the System Function Section. Eq. (2) shows that the 
transverse dipole field at the detector drops off as the cube of the 
distance from the voxel of tissue. By applying a voxel magnetizing 

10 field that increases as the cube of the distance from the detector to 
the voxel, the field at the detector becomes independent of the 
distance from the voxel. In this case, the field is given by 

B^W+fJ* 1 (26 ) 
where «, and B 0 arc constants and y ft is the distance from the detector 

15 to the voxel. Thus, a magnetizing field that increases as the distance 
cubed from the detector "levels" the magnitude of signals from voxels 
at different locations relative to the detectox. 



Detector Array 

20 

The ReMSI imager comprises a detector array of multiple 
detector elements which arc arranged in a plane. The array is a two 
dimensional detector array which is translated over the third 
dimension during the scan, or it is a three dimensional detector array. 

25 The individual detectors of the array respond 10 a single component 
of the secondary magnetic field which is produced by the magnetized 
object including tissue where the component of the field to which the 
detector is responsive determines the geometric system function 
which is used in the reconstruction algorithm discussed in the 

30 Reconstruction Algorithm Section. The detectors ideally point sample 
the secondary magnetic field at the Nyquist rate over the spatial 
dimensions which uniquely determine the susceptibility map which is 
reconstructed from the measurements. The magnetic field detectors 
include superconducting quantum interference devices (SQUlDs) 

35 wherein the SQUIDs are cooled below the critical temperature after 
the primary field is turned on and stabilized [6}. Thus, they link the 
flux of the primary field and are responsive to the secondary field 
only. Micromagnetic field sensors that are used to detect the primary 
field in the absence of the object to be scanned and to measure the 

40 secondary field have been developed that are based on galvanometric 
effects due to the Lorcntz force on charge carriers. In specific device 
configurations and operating conditions, the various galvanomagnetic 
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effects (Hall voltage, Lorcntz deflection, magnetorcsistive. and 
magnetoconcentration) emerge. Semiconductor magnetic field sensors 
include MAGFETs, magnetotransitors (MT), Van der Pauw devices, 
integrated bulk Hall devices including the vertical MT (VMT). and the 
5 lateral MT (LMT). silicon on saphire (SOS) and CMOS magnetodiodes. 
the magnetounijunction transistor (MUJT), and the carrier domain 
magnetometer (CDM), magnetic avalanche transistors (MAT), 
optoelectronic magnetic field sensors, and magnetoresistive magnetic 
field sensors. In the case of NMR measurement of the secondary field 
10 (and/or the primary field), the detector array comprises RF antennas 
described in the NMR Primary Magnet, Gradient Magnets, RF 
Generator, RF Transmitter, and RF Receiver Section. 

NMR Primary Magnet, Gradient Magnets, RF Generator. RF 
15 Transmitter, and RF Receiver 

The NMR embodiment of the RcMSI apparatus used to generate, 
and measure the secondary field and reconstruct the magnetic 
susceptibility map is shown in FIGURE &, and the corresponding 

20 coordinate system is shown in FIGURE 10. The apparatus comprises 
I.) a magnet including a superconducting magnet to magnetize a 
volume of an object or tissue to be imaged, 2.) a means 
(magnetometer) to determine the primary or magnetizing flux over 
the image space in the absence of the object to be imaged. In one 

25 embodiment, NMR of a proton containing homogeneous phantom of 
known susceptibility such as water is determined on a point by point 
basis to map the primary field, or the magnetic flux at multiple points 
is obtained simultaneously using the reconstruction algorithm 
described herein via determining the magnetic moment map over aii 

30 or a portion of the image space, and then determining the field using 
the known magnetic susceptibility of the phantom, 3.) a radio 
frequency (RF) generator and transmitter including an antennae such 
as saddle coils to excite the protons of the magnetized volume. 4.) a 
means including an antennae coil to sample the dipole component (z- 

35 component in terms of Eq. (2)) of the RF secondary magnetic field at 
the Nyqutst rate in time over the proton free induction decay, 5.) a 
detector array of elements of ihc means to sample the lime signals 
including an antennae array which is selectively responsive to the 
dipoic component (z~component of Eq. (2)) of the RF magnetic field of 

40 the magnetic moments of the protons which are aligned along the 

transverse axis to ideally point sample the secondary magnetic field 
at the Nyquisi rate over the spatiaS dimensions which uniquely 
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determine the susceptibility map which is reconstructed from the 
measurements. 6.) an analogue to digital convener to digitize the RF 
signals, 7.) a time Fourier transform processor to convert the signal at 
each detector over time into its spectrum, 8.) a reconstruction 
algorithm processor including a Fourier transform processor to 
convert frequencies into magnetic moments, map the magnetic 
moments to the image coordinates in parallel, convert each magnetic 
moment map into a magnetic susceptibility map, and superimpose the 
magnetic susceptibility maps of all the magnetic moments, and 9.) 
image processors and a display such that the ReMSI image can be 
rotated in space to be displayed from any perspective as a three 
dimensional or two dimensional (tomographic) image. 

The primary magnet is that typical of a magnetic resonance 
imaging apparatus. The primary magnetizing field is in the z- 
15 direction as shown in FIGURE 10. usually coaxial with the patient s 
body. However, in an another embodiment ihc primary maenetizing 
field is in the direction which is perpendicular to the patient's body. 
An open primary magnet design is preferable such as a Helmholtz coil 
design to accommodate ihe three dimensional RF antennae array. A 
20 further field is applied io have a gradient in the z-direction (FIGURE 
10). This provides a unique total value field in a chosen volume of 
the patient including a cross-sectional slice of the patient. A rotating 
RF field H„ of frequency chosen to cause resonance in the slice is then 
applied. Thus, only the molecules in the slice resonate. The 
25 resonance signal from the slice can then be delected. In another 
embodiment, the NMR (Larmor) frequency from each voxel is 
essentially determined by the magnitude of the magnetizing flux at 
the voxel position. Thus, the magnetizing flux provides a carrier 
frequency. The primary magnetic field has a sharp gradient at the 
30 edges of the image space wherein the Larmor frequencies outside of a 
selected range arc rejected as arising from voxels outside of the image 
volume. 

The ReMSI imager comprises a RF generator 14 and 18 and RF 
transmitter coils 9 shown in FIGURE 8 which provide the rotating II, 

35 (RF) field that is resonant with the protons of the magnetized volume 
and is oriented perpendicularly to the direction of the magnetizing 
field as shown in FIGURE 11. When the precessing nuclei are 
subjected io the additional rotating (RF) field H,. which is 
synchronous with the precession, their magnetic' moments and thus M 

40 prccesses about H, and rotate away from the primary field II 0 by an 
ancle 6 in a coordinate frame which rotates at the Larmor frequency 
[3}. The precession about H, continues as long as H, exists. The final 
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value of $ then depends on the strength of II,, which determines the 
precession rate, and the time for which it is turned on. In an 
embodiment, the strength and duration of H, is such that 0 = 90° such 
that the dipoie is oriented in the x R y R -planc of FIGURE II. A rotating 
5 coordinate frame is traditionally used to explain the physjes of NMR 
[3]. Thus, is terms of the traditional NMR coordinate designation as 
described by Patz [3] and Hounsfield [7] shown in FIGURE 11. the x R - 
axis and y R -axis rotate about the primary field H 0 oriented paral!ci\o 
the z R -axis at the Larmor frequency relative to the stationary NMR 

10 coordinate system (x, y, z) shown in FIGURE 10. Thus, the additional 
rotating (RF) field H, and the transverse RF magnetic field are 
stationary in the rotating NMR coordinate system, but rotate at the 
Larmor frequency in the stationary NMR coordinate system. Both the 
rotating and stationary coordinate systems are shown in FIGURE II. 

15 In an embodiment, the rotating H, (RF) field is along the y R -axis and 
rotates the magnetization vector by the angle 0 = 90° to into the x R y R ^ 
plane. In terms of the reconstruction by the present invention, the 
system function of the field corresponding to the rotated 
magnetization vector is equivalent to that of the detection of the z- 

20 component of a z-orientcd dipoie. The y~axis is the unique axis of the 
NMR system shown in FIGURE 8 and the stationary NMR coordinates 
shown in FIGURES 10 and II that corresponds to the unique axis of 
the Fourier Transform Reconstruction Algorithm Section, the z-axis 
shown in FIGURE 9, and the magnetization axis is the z-axis. The 

25 magnetization axis of the NMR rotating coordinates [3,7] is the z R -axis 
shown in FIGURE II. The measured transverse RF magnetic dipoie 
oriented in the xy-plane is periodically parallel to the y-axis and 
rotates about the primary field H, oriented parallel to the z-axis at 
the Larmor frequency. 

30 The RF transmitter includes saddle coils. The RF receiver is a 

three dimensional array, or two arrays 180° from each other, or two 
orthogonal pairs of arrays, wherein each member of a pair is 180° 
from each other. The detector ideally point samples the RF field at 
the Nyquist rate described in the Reconstruction Algorithm Section. 

35 In one embodiment shown in FIGURE 11, each detector 20a is a coil 
antennae perpendicular to the y-axis and selective to the y- 
componcnt of the RF secondary magnetic field arising from each voxel 
14a magnetized in the z-axis direction. This two dimensional array is 
translated along the y-axis during a scan where readings of the 

40 secondary magnetic field are obtained as a function of the translation. 
In another embodiment shown in FIGURE 12, the array 401 is three 
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dimensional comprising multiple parallel two dimensional arrays 402 
wherein each two dimensional array has a pluralily of antennae coils 
403 which detect the RF field from each magnetized voxel 404. 
In another embodiment, the rotating if, (RF) field shown in 
5 FIGURE II is along the x R -axis and rotates the magnetization vector 
by the angle <p = 9Q° to project into the x ft y R -plane. In terms of the 
reconstruction by the present invention, the system function of the 
field corresponding to the rotated magnetization vector is equivalent 
to that of the detection of the z-component of a z-oriented dtpole. The 

10 x-axis is the unique axis of the NMR system shown in FIGURE 8 and 
the stationary NMR coordinates shown in FIGURES 10 and 1! that 
corresponds to the unique axis of the Fourier Transform 
Reconstruction Algorithm Section, the z-axis shown in FIGURE 9, and 
the magnetization axis is the z-axis. The magnetization axis of the 

15 NMR rotating coordinates [3,7) is the z R -axis shown in FIGURE II. The 
detector array receives similar signals as those produced by the 
transmitter, and both are perpendicular to the x-axis. Each detector 
is a coil antennae perpendicular to the x axis and selective to the x~ 
component of the RF secondary magnetic field arising from voxels 

20 magnetized in the z-axis direction. This two dimensional array is 
translated along the x-axis during a scan where readings of the 
secondary magnetic field are obtained as a function of the translation. 
In another embodiment, the array is three dimensional comprising 
multiple parallel two dimensional arrays 402 as shown in FIGURE \2 

25 with the exception that the array 401 is perpendicular to the x-axis. 
The measured transverse RF magnetic dipole oriented in the xy-planc 
is periodically parallel to the x-axis and rotates about the primary 
field H 0 oriented parallel to the z-axis at the Larmor frequency. A 
linear combination of the cases of rotation of the RF magnetization 

30 along the x K ~axis and the y R -axis is within the scope of the present 
invention and may be adopted in a manner straightforward to those 
skilled in the NMR art to apply the fields and detectors described for 
this invention. 

The method described herein may be performed on a suitable 
35 NMR examining apparatus such as that shown in simplified form in 
FIGURE 8. and the corresponding coordinate system is shown in 
FIGURE 10. Illustrated schematically only are coils 6. which provide 
Z? 0 , the steady primary field 7, which provide G x > the field gradient in 
the x-axis direction as shown, 8 which provide the C r , (he field 
40 gradient in the y-axis direction as shown, 9 which provide the RF 
field, and 10. which provide G : , the field gradient in the z-axis 
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direction as shown. The coils arc driven by B 0 , G,, G v , RF, and C. 
drive amplifiers II, 12, 13, 14, and 15 respectively, controlled by B Qt 
G^, RF. and G control circuits 16, 17, 18, and 19, respectively. These 
circuits can take suitable forms which will be well known to those 
5 with experience of NMR equipment and other apparatus using coil 
induced magnetic fields. The circuits are controlled by a central 
processing and control unit 20 to achieve the desired primary field, 
field gradients, and RF field to rotate the magnetization vector such 
that it is perpendicular to coil 9. 

10 The RF coils are two saddle shaped coils 9 which are driven in 

parallel to provide the rotating RF field. The FID signals sensed are 
received in this example by the three dimensional array of RF coils 30 
and arc amplified by an RF amplifier 21 before being applied to signal 
handling circuits 22. The three dimensional detector array 30 is 

15 shown in more detail in FIGURE 12 as the three dimensional detector 
array 401. The circuits 22 arc arranged to make any appropriate 
calibrations and corrections but essentially transmit the signals to the 
processing circuits to provide the required representation of the 
examined volume. These circuits can conveniently be combined with 

20 the circuits which control the primary field, field gradients, and RF 

field and thus are included in the circuits indicated at 20, The picture 
thus obtained is viewed on a display 23, such as a television monitor, 
and this may include inputs and other peripherals 24 for the 
provision of commands and instructions to the machine, or other 

25 forms of output. 

The apparatus also includes field measurement and error signal 
circuits 25 which receive signals via amplifiers 26 from field probes 
X tt X }y Y %9 and Y : shown. 

The patient 27 is inserted in the tubular former of C, and G, 
30 coils 7, 8 and is supported there by a suitable couch or other 

supporting means. Such supports may be readily provided in any 
suitable form. 

The coils 7, 8 are two sets of coils axially displaced, each set 
comprising two pairs of saddle coils the pair 7 being at 90° to the pair 

35 8. These coils are themselves inserted into the central aperture in B 0 
coils 6 which in an embodiment are wound in four parts connected in 
series to provide an approximately circular configuration which is 
well known to be desirable for production of a uniform field. Further 
details of the coil winding will not be given since suitable coils can 

40 readily be devised, by those with the appropriate skills, to provide 
the fields required. 
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The appropriate stores provide the amplitude and duration 
signals which are converted to analogue form in digital to analogue 
converters (DACs) and applied to respective coil drive circuits x, y T z, 
RF. The respective drive circuits, which can take any form well 
5 known for driving field coils, provide the specified currents to the 

appropriate coil for the specified duration. The apparatus and circuits 
described so far may be adopted to provide different gradients and 
RF fields, by appropriately adjusting the stored sequences and profile 
data. Similarly other known NMR apparatus which are capable of 

10 applying a steady magnetic field, a pulsed RF field, and C„ G r . and C. 
field gradients to a body, may be adopted in a manner 
Straightforward to those skilled in the NMR art to apply the fields 
described for this invention. 

With the basic signal handling system of the present invention. 

15 the FID signals from the signal sensing coils of the detector array 30* 
(shown in more detail in FIGURE 12 as coils 403 of the three 
dimensional detector array 401) are amplified in an RF amplifier and 
applied via an analogue to digital converter (ADC) to a store such as a 
random access memory (RAM). The data is then processed according 

20 to the procedure given in the Reconstruction Algorithm Section. 

Reconstruction Algorithm 

In the reconstruction process described herein, the secondary 

25 field may be in Ihe same or in a transverse orientation relative to the 
primary field. In each case, the orientation of the secondary dipolc 
field and the measured secondary dipole field component are 
according to Eq. (2). The relationship of the coordinate systems of the 
present invention are shown in FIGURES 9, 10, and 11. The primary 

30 and secondary fields are parallel and stationary in FIGURE 9 versus 
transverse in the rotating and stationary NMR coordinates shown in 
FIGURES 10, II. and 12. The relationship of the coordinate system of 
the reconstruction method of the NMR system used herein to that of 
the coordinate system of Eq. (2) shown in FIGURE 9 is according to 

35 FIGURE II. The magnetization axis of Eq (2) is the z-axis; whereas, 

the magnetization axis of the stationary NMR coordinates is the z-axis 
shown in FIGURE 10, and the magnetization axis of the NMR rotating 
coordinates {3/7) is the z R -axis shown in FIGURE II. Regarding Eq. (2) 
and the reconstruction algorithm in the NMR case, the z-component of 

40 Eq. (2) is substituted with the stationary NMR y-component, the y- 
component is substituted with the stationary NMR z-component, and 
the x-component is substituted with the stationary NMR x- 
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component. For example, in the NMR case shown in FIGURES 8, 10. 
IL and 12, Step 7 of the Reconstruction Algorithm of the Fourier 
Transform Reconstruction Algorithm Section is equivalent to the 
general case except for the coordinate designations and the 
5 relationship between the magnetic moment of a static secondary field 
and the bulk magnetization M of the RF field wherein m.. the 
magnetic moment along the z-axis. of Eq. (2) corresponds to the bulk 
magnetization M of each voxel: 



10 7 > D,v * d e each element of Matrix C by to correct for the 

restriction that the sample space is defined as y greater than 
zero, y>0. This operation creates Matrix D which is the bulk 
magnetization M map. 

= D(x % y\z) 

15 ln one embodiment, the voxels with a finite bulk magnetization 

M above a certain threshold or at an edge in the map arc 
assigned the magnetic susceptibility (magnetic moment) 
determined from the shift of the Larmor frequency due to the 
presence of the voxel(s) in the uniform (nonuniform) 

20 magnetizing field. The other voxels are assigned a zero value. 

This procedure is repeated for all Larmor frequencies, In the 
limit with sufficient Larmor frequency resolution, each volume 
element is reconstructed independently in parallel with all 
other volume elements such that the scan time is no greater 

25 than the nuclear free induction decay (FID) time. In the case 

that the primary field is nonuniform, the magnetic moment map 
is determined, then the magnetic susceptibility map of each 
magnetic moment is given by dividing the magnetic moment 
map by the primary magnetic field map on a voxel by voxel 

30 basis and by subtracting one from each term (Eq. (23)). The 

total magnetic susceptibility map is the superposition of the 
separate maps for each magnetic moment corresponding to a 
unique Larmor frequency which is plotted and displayed. 

35 In summary, the magnetic moment is known from the Larmor 

frequency. The NMR method then determines the voxel position(s) 
using the measurements of the spatial variations of RF field wherein 
each set of measurements over the detector array corresponding to a 
given frequency can be normalized to some convenient magnitude. 
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The bulk magnetization M map is transformed into a magnetic 
susceptibility map. 

The reconstruction algorithm can be a reiterative, a matrix 
inversion, or a Fourier Transform algorithm. For all reconstruction 
5 algorithms, the volume to be imaged is divided into volume elements 
called voxels, and the magnetized voxel with magnetic moment z& is 
modeled as a magnetic dipole where x »s the magnetic susceptibility 
of the voxel and B is the magnetizing flux at the position of the dipole. 
The reconstruction may be multiplexed by using the NMR signals as 
10 the input. The transverse RF magnetic field of each voxel is that of a 
dipole, the maximum amplitude is given by Eq. (2) wherein the 
Larmor frequency of each voxel is shifted due to its magnetic 
susceptibility, and m : , the magnetic moment along the z-axis, of Eq. 
(2) corresponds to the bulk magnetization M of each voxel. The NMR 
15 frequency is determined by the magnetizing flux which includes a 

contribution from the magnetic susceptibility of the voxel. Each voxel 
or voxels which have a unique local magnetic flux will have a unique 
a NMR frequency. The spatial variation of the intensity of the 
transverse RF signal of this unique frequency is the input for the 

20 reconstruction algorithm. This procedure is performed either in 
parallel or series for each data set of signals of each unique 
frequency. The magnetic susceptibility image is the superposition of 
the independent images, each of which corresponds to a unique 
frequency. The general steps are: |.) The magnetic susceptibility 

25 (magnetic moment) of each voxel(s) is determined from the shift of 
the Larmor frequency due to the presence of the voxei(s) in the 
uniform (nonuniform) magnetizing field. 2.) The measurements of the 
spatial variations of the transverse RF field of a given frequency is 
used to determine the coordinate location of each voxel(s). Thus, each 

30 matrix comprises the intensity variation over the sample space of the 
RF field of the bulk magnetization M of each voxel having the 
corresponding magnetic moment. 3.) In the case that the magnetizing 
field is uniform, the magnetic susceptibility map is determined 
directly, but in the case thai the primary field is nonuniform, the 

35 magnetic moment map is determined, then the magnetic 

susceptibility map of each magnetic moment is given by dividing the 
magnetic moment map by (he primary magnetic field map on a voxel 
by voxel basis and by subtracting one from each term (Eq. (23)). 4.) 
The total magnetic susceptibility map is the superposition of the 

40 separate maps for each magnetic moment corresponding to a unique 
Larmor frequency which is plotted and displayed. 
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A matrix inversion reconstruction algorithm is to use the 
geometric system function to determine the spatial intensity variation 
of the transverse RF field over the detector array, to invert the 
corresponding matrix, and to multiply the signal over the detector 
5 array by the inverted matrix to give the voxel sources. For example, 
a matrix inversion reconstruction algorithm is to determine a 
coefficient for each voxel mathematically (Eq. (30)) or by calibration 
which when multiplied by the bulk magnetization M of each voxel is 
that voxel's contribution to the signal at a given detector at a given 

10 Larmor frequency. This is repeated for every detector, and those 

coefficients are used to determine a matrix which, when multiplied by 
a column vector of the bulk magnetization M values of the voxels, 
gives the voltage signals at the detectors. This matrix is inverted and 
stored in memory. Voltages as a function of time are recorded over 

15 the detector array. The signal as a function of time is Fourier 

transformed to give the spectrum. The intensity variations over the 
detector array at each Larmor frequency are multiplied by the 
inverse matrix, to generate the bulk magnetization M map. In one 
embodunent, the voxels with a finite bulk magnetization M above a 

20 certain threshold or at an edge in the map are assigned the magnetic 
susceptibility (magnetic moment) determined from the shift of the 
Larmor frequency due to the presence of the voxel(s) in the uniform 
(nonuniform) magnetizing field. The other voxels are assigned a zero 
value. This procedure is repeated for all Larmor frequencies. In the 

25 limit with sufficient Larmor frequency resolution, each volume 
element is reconstructed independently in parallel with all other 
volume elements such that the scan time is no greater than the 
nuclear free induction decay (FID) time. In the case that the primary 
field is nonuniform, the magnetic moment map is determined, then 

30 the magnetic susceptibility map of each magnetic moment is given by 
dividing the magnetic moment map by the primary magnetic field 
map on a voxel by voxel basis and by subtracting one from each term 
(Eq. (23)). The total magnetic susceptibility map is the superposition 
of the separate maps for each magnetic moment corresponding to a 

35 unique Larmor frequency which is plotted and displayed. 

A reiterative reconstruction algorithm is to use the geometric 
system function to determine the system of linear equations which 
gives the spatial intensity variation of the transverse RF field over the 
detector array for each Larmor frequency. Voltages as a function of 

40 time are recorded over the detector array. The signal as a function of 
time is Fourier transformed to give the spectrum. The sources of the 
RF field are the voxels of the image space each having a magnetic 
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susceptibility x and a magnetic moment contribution // 0 jlf (Eq. (22)) 
which determines the Larmor frequency, and a bulk magnetization M 
at a given location which determines the intensity and spatial 
variation of the RF field over the sample space. The system of linear 
5 equations relates the voltage at each frequency at each detector to the 
bulk magnetization M of each voxel which contributes signal to the 
detector at a given Larmor frequency due to its magnetic 
susceptibility x and magnetic moment contribution /i^H (Eq. (22)). 
For each Laimor frequency, each of these equations is the sum over 

10 all voxels of the bulk magnetization M value of each voxel times its 
weighting coefficient for a given detector which is set equal to the 
voltage at the given detector. The coefficients are determined 
mathematically (Eqs. (30), (21), (22), and (23)), or they are 
determined by calibration. The bulk magnetization M for each voxel 

15 is estimated, and the signals at each detector are calculated using 
these values in the linear equations. The calculated values are 
compared to the scanned values and a correction is made to M of each 
voxel. This gives rise to a second, or recomputed, estimate for M of 
each voxel at each Larmor frequency. The signal value from this 

20 second estimate is computed and corrections are made as previously 
described. This is repeated until the correction for each reiteration 
approaches a predefined limit which serves to indicate that the 
reconstruction is within reasonable limits of error. The reiterative 
procedure is repeated for all Larmor frequencies. The total magnetic 

25 susceptibility map is the superposition of the separate maps for each 
magnetic moment corresponding to a unique Larmor frequency which 
is plotted and displayed. 

The general process of reconstruction by reiteration is shown 
according to the steps of FIGURE 13 (and is implemented in processor 

30 20 in FIGURE 8). The image displayed according to the process 200 is 
directly related to the bulk magnetization M of voxel sections of the 
object examined, the image is merely a mapping of the magnetic 
susceptibility in three dimensions which determines the Larmor 
frequency. Accordingly, signals produced by the RF magnetic sensors 

35 110, in terms of volts as a function of Larmor frequency, are a direct 
result of the bulk magnetization M of the voxel elements. Therefore, 
a reference voltage is generated at 210 from which the actual or 
measured sensor voltages is subtracted at 220. The reference 
voltages are modeled by assuming a signal contribution from each 

40 voxel element to each sensor. Therefore, the signals are summed 

separately for each sensor at each frequency according to a weighted 
contribution of the voxel elements over the x, y, and z-axes in terms 
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of the coordinates of Eq. (2) and FIGURE 9 according to a model of the 
tissue to be examined. The resulting modeled or calculated voltage 
signals are compared at step 220, providing a difference or A signal, 
weighted at step 230 to produce a weighted difference signal, which 
5 is then added to (he previously estimated bulk magnetization M 
value for each voxel element at step 240. The resulting level, 
available in three dimensions corresponding to the axes x. y, and z, is 
selectively displayed on the display at step 250. Having adjusted the 
estimated bulk magnetization M for each voxel, the calculated bulk 

10 magnetization M is recalculated at step 260, the resulting estimated 
sensor voltage is then compared to the actual sensor voltage at step 
220, the process 200 being repeated until the difference is reduced to 
a tolerable value. In one embodiment, as the final step or as input to 
step 220. the voxels with a finite bulk magnetization M above a 

15 certain threshold or at an edge in the map are assigned the magnetic 
susceptibility (magnetic moment) determined from the shift of the 
Larmor frequency due to the presence of the voxei(s) in the uniform 
(nonuniform) magnetizing field. The other voxels are assigned a zero 
value. This procedure is repeated for all Larmor frequencies. In the 

20 limit with sufficient Larmor frequency resolution, each volume 
element is reconstructed independently in parallel with all other 
volume elements such that the scan time is no greater than the 
nuclear free induction decay (FID) time. In the case that the primary 
field is nonuniform, the magnetic moment map is determined, then 

25 the magnetic susceptibility map of each magnetic moment is given by 
dividing the magnetic moment map by the primary magnetic field 
map on a voxel by voxel basis and by subtracting one from each term 
(Eq. (23)). The total magnetic susceptibility map is the superposition 
of the separate maps for each magnetic moment corresponding to a 

30 unique Larmor frequency which is plotted and displayed. 

In terms of the coordinates of Eq. (2) and FIGURE 9. the 
reconstruction algorithm using Fourier Transforms (FFT) involves 
exploiting the FFT to solve Eq. (30) given below. The FFT algorithm is 
given below and is followed by a discussion of the resolution of the 

35 magnetic susceptibility map in the Nyquist Theorem with the 

Determination of the Spatial Resolution Section. For the case that 
follows, data is acquired in the x, y t and z-directions, but in general, 
data is acquired over the dimensions which uniquely determine the 
magnetic susceptibility map. Also, the present analysis is for 

40 measuring the z-component of the magnetic field of a dipole oriented 
in the z-direction; however, the analysis applies to the other two 
orthogonal components where the geometric system function for the 
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z-cornponent of the z-oriented dipole is replaced by the geometric 
system function for the x or y-componems of the magnetic field 
produced by the dipole where the geometric system function is 
defined below as the impulse response of the detector to the given 
5 component of the Held of a dipole of given orientation. The sample 
space, or space over which the secondary RF field is measured, is 
defined in the present example as the three-dimensional space 
comprising the entire xy-plane and the positive z-axis, as shown in 
FIGURE 9. Other sample spaces are valid and each requires special 

10 consideration during the reconstruction as described below. The 

discrete voltages recorded from an infinite detector array in the xy- 
planc which is translated to infinity along the z-axis starting from the 
origin where the detector array is responsive to the z-cornponent of 
the secondary magnetic field is given by Eq. (2), where the voltage at 

15 any point in space produced by dipoles advanced in the z-direction 
and advanced or delayed in the x and y-directions relative to the 
origin is given by the following Eq. (29), where the voltage response is 
Q times the secondary magnetic flux strength in the case shown in 
FIGURE 9. 

20 The spatial variation of the normalized RF signal at each Larmor 

frequency is used in the NMR case shown in FIGURE 11. The 
relationship of the coordinate designations of Eq. (2) and FIGURE 9 
versus the NMR system and the relationship between the magnetic 
moment of a static secondary field and the bulk magnetization M of 

25 the RF field is given in the Reconstruction Algorithm Section and in 
the Fourier Transform Reconstruction Algorithm which follows. 

Fourier Transform Reconstruction Algorithm 

30 In terms of the coordinates of Eq. (2) and FIGURE 9, the- volume 

to be imaged is divided into volume elements called voxels and the 
magnetized voxel with magnetic moment density ^B is modeled as a 
magnetic dipole where x *s *he average volume magnetic 
susceptibility of the voxel and B is the magnetizing flux at the 

35 position of the dipole. Data comprising the z-cornponent of the 

magnetic field of a dipole oriented in the z-direction is acquired in the 
three dimensional sample space comprising the xy-plane and the 
positive z-axis. The magnetic moment, /n_, of each voxel within the 
phantom is given by the product of its volume magnetic 

40 susceptibility, x> &t magnetizing flux oriented along the z-axis, B ;> 
and the volume of the voxel, Av\ 

in. =^BAV (27) 
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The magnetic moment of each voxel is a magnetic dipole. And the 
phantom can be considered to be a three-dimensional array of 
magnetic dipoles. At any point extrinsic to the phantom, the z- 
component of the secondary flux t B, from any single voxel is 

B ^ m :TT~~* ^577 (28) 

(or + v +z ) 

where x, y, and z are the distances from the center of the voxel to the 
sampling point. It is shown in APPENDICES MV that no geometric 
distribution of magnetic dipoles can give rise to Eq. (2). Therefore, 
the flux of each magnetic dipole (voxel contribution) forms a basis set 
for the flux of the array of dipoles which comprise the magnetic 
susceptibility map of the phantom. 

Eq. (2) is a system function which gives the magnetic flux 
output in response to a magnetic dipole input at the origin. The 
phantom is an array of spatially advanced and delayed dipoles 
15 weighted according to the magnetic susceptibility of each voxel; this is 
the input function. The secondary flux is the superposition of 
spatially advanced and delayed flux, according to Eq. (2); this is the 
output function. Thus, the response of space to a magnetized 
phantom is given by the convolution of Eq. (2) with the scries of 
20 weighted, spatially advanced and delayed dipoles representing the 

susceptibility map of the phantom. The discrete voltages recorded by 
a detector array over the sample space comprising the xy-plane and 
the positive z-axis are 

i w i=c/f f f Aj , 4»^: l '- |y -f: l< -^l ( 29, 

»>*o . : -, ;/ >i ^- U/ 2> }(*-«,*,] +[y~n I ib } ) +[* + ",*>]] 

25 where the flux magnetizing each voxel is unity, the volume element is 
made unity, and the voltage response of each detector is C 0 times the 
secondary magnetic flux strength. The variables of Eq. (29) are 
defined as follows: 

30 the magnetic susceptibility of the voxel 

located at the position defined by the 
Dirac delta function, 

8{x-n { k t , y-n : k z . Z-n 5 k y ) 

35 dipole spacing in the x t y, and z- 

directions, respectively 
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/ Jt />, /j the dimensions of the image space, the 

total distance in the x, y. and z- 
directions, respectively, for which the 
magnetic susceptibility of the voxels is 
5 nonzero 

s y the detector spacing in the x, y, and z- 

directions, respectively 

10 The voltage signals recorded at the detector array over the 

sample space is 

V(m |V m 2i m,]s Z Z Z C A X " m ^y-^2 s i^ ~ m ) s >] 



(30) 



L I, Z n — 1 — ^>®x* t « .«A* n : k *- z - n M 

Eq. (30) can be represented symbolically as follows 

s = C 0 \gx\h®f\x»{z)] (31) 
15 where C 0 is the proportionality constant between the signal voltage 
and the output flux strength, s is the discrete function of voltage 
signals recorded of the secondary flux over the sample space, g is the 
secondary magnetic flux sampling function given as follows: 

S^Z Z Z ^-v P >-«^^ (32) 

20 In Eq, (31), h is the system function which is also defined as the 
geometric system function given as follows: 

7~~ i ; TTvT (33) 
(*- + y* -h*-) 

The system function is the impulse response of the detector array 
where the external magnetizing flux is set equal to unity. [If the 

25 magnetizing field is nonuniform, then the magnetic susceptibility 
(magnetic moment) values determined from this Fourier transform 
algorithm (Eq. (62)) must be divided by the strength of the magnetic 
flux magnetizing the corresponding voxel on a voxel by voxel basis. 
This procedure is described in detail in the Altering the Dynamic 

30 Range Provided by the System Function Section.] In Eq. (31), / is the 
magnetic susceptibility function given as 

' = Z Z Z x^„A*-*%kvy-*fo*i-*&\ ( 34 > 
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In Eq. (31), u(z) is the unitary z function which is one for positive z 
and zero otherwise. The function g discretizes the continuous voltage 
function V given by Eq. (29) which is the function h convolved with 
the function / then multiplied by the function n(z). The discrete 
voltages recorded over the sample space are used in a Fourier 
transform algorithm derived herein to reconstruct the magnetic 
susceptibility image. Consider, the function s of Eq. (30) which is 
given as follows: 

s~{h®f)xu{i) (35) 
Eq. (35) is equivalent to the function h convolved with the function / 
then multiplied by the function u(z). The function S which is the 
Fourier transform of s is given as follows- 

S = (HxF)*U(k : ) (36) 
The function 5 is equivalent to the resultant function of the function 
H. the Fourier transform of the system function- ft, multiplied by the 
function F, the Fourier transform of the magnetic susceptibility 
function--/, convolved with the function */(*.), the Fourier transform 
of the unitary z function- u{i). The Fourier transform of the magnetic 
susceptibility function -/ (Eq. (34)) 
•«./» .<-/»( .r,/ji 

t0 /=Z Z X X,, „..., S[r y z-n,kA (37) 



15 



30 



is 



f =Z I II e-*-— I (38 ) 

where jr. = /!,*,; y a = Z(i = The Fourier transform of 

u{z)=\forz>Oandu(z)=\forz<0 is (8J 

25 (39) 

The Fourier transform of the system function (See APPENDrX V) 

(, l + y ! + ^f 2 "(p> + ^ f , < 4 <» 



IS 



4 a*; 



where 
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* r = 2;rf, = 2;r- 
y 

I 



(43) 



k. = 2nf.=2x- 
z 

The resultant function of the product of the functions H (Eq. (41)) and 
F (Eq. (38)) then convolved with the function U(k ) (Eq. (39)) is 
5 = HF ® (/[*.} = 

2nkl -'• 



+ T®F7f£ X Ii. '-i 



(44) 



5 The result is given as follows, and the derivation appears in 
APPENDIX VI. 

*=^X Z Iw 



Substitution of the Fourier transform of the system function, H (Eq. 
10 (41)), the Fourier transform of the magnetic susceptibility function. F 
(Eq. (38)), and factoring out * p l in the second term gives 

-', 



(45) 



5 = 



2nk\ 2«*j 



EES 



W : /2 



+ X Z Z ^w*-"*-* 

:.«0 v, »-/./» «.»-i,/2i- 



i+y- 



Multiplication by the complex conjugates gives 



(46) 
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7 Lu A—t A-n, 

^1.-0 >.— UA X.*-U2l 



Factoring out the common terms gives 

K k l _~b •<:/* ♦'»/* 



l 


ht] 


P . 


ht] 



1 









(47) 



4*«*« **»>'- **;*•{ 



(48) 



Combining exponential terms gives 
5 



The function S divided by the function // is 

" y.*-/./2 [ Ar p t 



(49) 



«.<0 

The inverse Fourier iransform of the function S divided by the 
10 function H is given as follows, where the symbol F % {Q) is defined as 
the inverse Fourier transform of the function Q . 

-y y y y ,->K'.**.>.*^-i 

Q W l..* t 



(50) 



O ^ 



c.<0 



(51) 



The inverse Fourier transform of the function F is 

= IXz w «|»-*,)i-y..:-;] (52) 

y, j. 

15 Specifically, the inverse Fourier transform of the function F in (erms 
of the dipole spacing is 
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Z 5j Z £«i.'»j.n J < 

»\*0 J.— i,=-<,/2 



(53) 



= Z Z Z ^,,^*-^y-^-"3*) 

Consider the general Fourier transform pair in cylindrical coordinates: 

*M.z)»GM.fl>) (54) 
Under circular symmetry, that is when £ is independent of 9 (and 
5 hence G independent of #) the inverse Fourier transform is [9] 



10 



15 



h(r.z) = 2/rJ jHi^co^ilKs^e^sdsdo) 
o — 

The inverse Fourier transform of Eq. (51): 

n Z-r Z-r Z-f *<i|.«».*, C 

given by Fq. (53) and Eq. (55) is 



j 



(55) 



*„<0 



(56) 



2*j j ^'''^(vlVV'*'^* 
+/2* j J e'^jjlyi)*,*/* : dk. 



— o 



The solution of Eq. (57) appears as follows, and the derivation appears 
in APPENDIX VII. 

t 1 I *. "4 e"^ 1 + " ♦ j _ ,-M,l] 

l^'.-O r, a -l-/l jr.^/i [ V ' 



(57) 



=f Z Z Z x^A*-"^-'^-^]® 



Inverse Transform I 



1ST 



(58) 
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Thus, the inverse Fourier transform of — (Eq. (51)) is given by Eq. 
(58). 



-fx 



W./J4 ♦/,/}! 



f~~i JOT 



The convolution replaces each coordinate with a spatially shifted 
coordinate. 



10 The result of the evaluation of Eq. (60) at each coordinate x^y^ is 
Solving for the magnetic susceptibility of each dipole gives 

-[11 



It 

2zl 



(59) 



(60) 



(61) 



(62) 



The solution of the magnetic susceptibility of each dipole of the 
15 magnetic susceptibility function follows from Eq. (62), Discrete values 
of the voltages produced at the detector array due to the secondary 
magnetic field are recorded during a scan which represent discrete 
values of function 5 (Eqs. (30-31)); thus, in the reconstruction 
algorithm that follows, discrete Fourier and Inverse Fourier 
20 transforms replace the corresponding continuous functions of Eq. (62) 
of the previous analysis- 
Discrete values of H of Eqs. (41-43), the Fourier transform of 
the system function, replace the values of the continuous function. 
Furthermore, each sample voltage of an actual scan is not truly a- 
25 point sample, but is equivalent to a sample and hold which is 

obtained by inverting the grid antennae matrices or which is read 
directly from a micro-antennae as described in the Finite Detector 
Length Section. The spectrum of a function discretely recorded as 
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values, each of which is equivalent to a sample and hold, can be 
converted to the spectrum of the function discretely recorded as point 
samples by dividing the former spectrum by an appropriate sine 
function. This operation is performed and is described in detail in the 
Finite Detector Length Section. From these calculated point samples, 
the magnetic susceptibility function is obtained following the 
operations of Eq. (62) as given below. 

Reconstruction Mgft"th m 
1 ) Record the voltages over the sample space. 



2) Invert the detector grid matrices defined by the 

noncommon areas of the overlapping elements of the 
15 detector array described in the Finite Detector Length 

Section to obtain the sample and hold voltages which form 
Matrix A (if micro-antennas are used, form Matrix A of 
the recorded voltages directly). 

20 3 ) Three-dimensionally Fourier transform Matrix A. using a 

discrete Fourier transform formula such as that which 
appears in McC Siebcrt [101 to form Matrix A* of elements 
at frequencies corresponding to the spatial recording 
interval in each direction. 

4) Multiply each element of Matrix A* by a value which is 
the inverse of the Fourier transform of a square wave 
evaluated at the same frequency as the element where 

30 *he square wave corresponds to a sample and hold 

operation performed on the continuous voltage function 
produced at the detector array by the secondary field. 
This multiplication forms matrix A** which is the discrete 
spectrum of the continuous voltage function discretely 

35 point sampled (See the Finite Detector Length Section for 

_ details of this operation). 

A * {k <- k > - k)x 7- fJi n 55 A * *<Mr. *. ) 

sinc{*, .*,.*.) 

5) Multiply each element of Matrix A** by the value which is 
40 * h e inverse (reciprocal) of the Fourier transform of the 
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20 



system function evaluated at the same frequency as the 
element to form Matrix B. 

A JUL)* = B(k k k) 

6) Inverse three-dimensionally Fourier transform Matrix B 
using ihe discrete inverse Fourier transform formula such 
as that which appears in McC Siebert [10} to form Matrix 
C whose elements correspond to the magnetic 
susceptibility of the dipoles at the points of the image 
space of spatial interval appropriate for the frequency 
spacing of points of Matrix B, 

7) Divide each element of Matrix C by — to correct for the 

restriction that the sample space is defined as z greater 
ihan zero (z>0). This operation creates Matrix D *which is 
the magnetic susceptibility map 
C(x.y.z) 

= D(x.y,z) 

(If the magnetizing field is not constant, then further 
divide each element by the value of the magnetizing field 
at the position of the corresponding magnetic 
susceptibility element.) 



7*) In the NMR case shown in FIGURES 8, 10, ll, and 12. Step 
25 7 of l nc Reconstruction Algorithm of the Fourier 

Transform Reconstruction Algorithm Section is equivalent 
to the general case except for the coordinate designations 
and the relationship between the magnetic moment of a 
slalic secondary field and the bulk magnetization M of the 
30 RF field: 

Divide each element of Matrix C by ~ to correct Tor the 

restriction thai the sample space is defined as y greater 
than zero, y>0. This operation creates Matrix D which is 
35 the bulk magnetization M map. 
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C{x,y f z) 

£ = 0(*,y.2) 

(In one embodiment, ihe voxels with a finiie bulk 
magnetization M above a certain threshold or at an edge 
in the map are assigned the magnetic susceptibility 
(magnetic moment) determined from the shift of the 
Larmor frequency due to the presence of the voxei(s) in 
the uniform (nonuniform) magnetizing field. The other 
voxels are assigned a zero value. This procedure is 
repeated for all Larmor frequencies. In the limit with 
sufficient Larmor frequency resolution, each volume 
element is reconstructed independently in parallel with all 
other volume elements such that the scan time is no 
greater than the nuclear free induction decay (FID) time. 
In the case that the primary field is nonuniform, the 
15 magnetic moment map is determined, then the magnetic 

susceptibility map of each magnetic moment is given by 
dividing the magnetic moment map by the primary 
magnetic field map on a voxel by voxel basis and by- 
subtracting one from each term (Eq. (23)). 

20 

8) Plot the magnetic susceptibility function with the same 

spatial interval in each direction as the sampling interval 
in each direction. 

25 8 *) In lh e NMR case shown in FIGURES 8. 10, IK and 12. Step 

8 of the Reconstruction Algorithm of the Fourier 
Transform Reconstruction Algorithm Section is equivalent 
to the general case except for the coordinate designations, 
the relationship between the magnetic moment of a static 

30 secondary field and the bulk magnetization M of the RF 

field, and the determination of the magnetic susceptibility 
of the mapped voxels by the Larmor frequency (Eqs. (20- 
23)) in the NMR case: 

Plot the magnetic susceptibility function 

35 corresponding to the Larmor frequency for the mapped 

voxels with the same spatial interval in each direction as 
the sampling interval in each direction over all Larmor 
frequencies (i.e. the total magnetic susceptibility map is 
the superposiiion of the separate maps for each magnetic 
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moment corresponding to a unique Larroor frequency 
which is plotted and displayed). 

(The above steps for each Larmor frequency relate generally to the 
program implementation shown in the listing of the Exemplary 
Reconstruction Program as follows. The above Steps 1 and 2 relate to 
(he Data Statements beginning at lines 50; and Step 3 relates to the X 
Z and Y FFT operations of lines 254, 455 and 972, respectively. Steps' 
4 and 5 are implemented by the processes of lines 2050, 2155 and 
2340; and Step 6 relates to the X, Y and Z inverse transform of lines 
3170, 3422 and 3580, respectively. Step 7 and 7* relates to the 
correction and normalization process of line 4677.) 

Exemplary Reconstruction Program 



10! ReMSI ALGORITHM 

25 OPEN "ReMSI PROTOTYPE I .LIS" FOR OUTPUT AS # 1 

26 C=.0S -• 

20 27 DIM X(9.9,9) 

28 FOR, W=l TO 9 STEP 1 

29 FOR T= I TO 9 STEP 1 

30 FOR S = l TO 9 STEP 1 

31 X(W,T,S>^0 
25 3 2 NEXT S 

3 3 NEXT T 

34 NEXT W 

35 X<S.5.6J=1 

36 DIM DII9.9) 

30 37 PRINT II, -DIPOLE PHANTOM " 

40 FOR 0=1 TO 9 STEP I 

41 FOR R=l TO 9 STEP 1 

4 2 FOR U=l TO 9 STEP 1 
43 LET H=X(U.R,0> 

35 44 LET DI(U.R)=H 

45 NEXT U 

46 NEXT H 

47 MAT PRINT A 1 . DI . 

48 NEXT Q 

40 50 DATA 0,0,0.0.0,0,0,0,0,0,0,0,0 

52 DATA 0.0. 0.0.0.-. 6.-. 9. -I. 1.-1.2. -1.8. S. 7. 3 5 4 - 4 - 5 - 7 - 6 7 
54 WTA 4.2.19.7,7.6.. 1...9.-.3... 4... 7 . ; 7 .ii 9 2S 6S J Ji . i ! 1 9 J 1 1 8 J ! 1 3! . 

45 i 6 o DATA '• 4 - 2 - 7 ' 10 - 7 - 12 - 4 - 5 -^^ 

58 DATA 

^ O ;^'^^ O ' 0 '°' 0 ' O - O - 0 ' O '0'0,0.0.0,0. 0.0. 0.0.0.0. 0.0.0.0.0,0.0.0 
60 DATA 0,0, -.5. -.8, -.8. -1.0,-1.2, - .8.6.9, 2^2-4-1 
.6.1.2.4.1.12.7,0, .6 
50 62 DATA -.6. -.3. -.5. -.7. -.3, 3. 7. 9. 05, 7. 6, 1.3 - 5 - 3 - 4 ~ 
.2,2.7.4.3. 10.6. .8 ' 
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64 DATA 1-7. -.4, -.3.. .5. -.4. -.1.1. 1.1. 1,1. 7, .2.- 

66 DATA 0,0. 0.0,0. 0.0. 0.0. 0,0. 0.0. 0.0. 0.0. 0.0, 0.0. 0,0. -.5. -.6, -.6, 

68 DATA .3, 2. 5, .l.-.l.-. 3, -.4. -.3, 1-3. 3. 9. 8. 1.1. 7. .9. -.3. -.3.- 4 - 6 
.8.1.5 " r 

? 2 Tl ^ 6S - 4 ' 7 ' 1 - 3 ^- 3 --' 2 ^ 2 - 1 -2^^.5.7.1..7,.7,-.3,>.2. .4,- 

72 DATA 1.1, 1.4. .4.- 

-3.0,0,0,0.0,0,0,0,0.0,0,0,0.0.0.0,0,0,0,0,0,0.0,0,0,0 

74 DATA 0,0, 0,0. 0,0. 0.0, 0,0, -.4, -.5,0 -1-11715 2-1-2 

-2,0.1.3.3.1 ' ' 

™ DATA 4. 1.3. 4. .9,0,-. 2, -.3, -.4, -.4,1.1. .6S, 2. 9. 1.1. -.2 - i - 
1. .2,1.2.1.6 

78 DATA 4. 9, .1.. 5, -.3, -.2, -.3, -.1.0, .3. 1.3. .6 1 - 
• 1.0.0.0.0.0,0,0,0,0.0,0.0 

80 DATA 0.0.0,0,0,0,0.0,0,0,0.0,0.0.0.0,0.0,0,0,0.0,0 0 - 3 - 
•2..1, .9, 1.6, 1.6 • . • . 

82 DATA .4.O.O.-.I. -.1. .2. 1.0,2.2.2.7.3.7, .5 O.-.2.- 2 - 2 - 
-1, .6, .25,1.6, .7 

84 DATA ~.2,-.l.-_l. .1..7, .9.2.8,0, . 3 . - . 2 , - . 1 , - . 2 . - . 1 , 0 . 2 3 3 0 - 
•1,0,0,0 

85 DATA 0.0.0,0,0.0,0,0.0,0.0.0.0.0,0 

86 DATA G. 0.0. 0.0. 0,0.0. 0,0, 0.0, 0,0, 0.0,0 0 - 2 ~ 
.1., i.l.0. 1.4. 1.3.. 2. 0.6 

88 DATA -.1,0, .2, .7,1.5.2.0,2.9. .3.0,-. l,-.l, .1,0 7 25 9 4 - 

.1,0,0..!. .5 

90 DATA .5.1.3.0. .2. -.1.0. -.1.0. .1, .2,. 1. .2.0. - 
.1, 0,0,0.0. 0,0, 0,0. 0.0. 0.0, 0.0 

92 DATA 0.0,0.0,0,0.0.0.0.0,0,0.0.0,0,0.0 0 0 0 0 0- 
-1,0, .1, .8, .9.1.0, .1.0.0 

94 DATA 0,0. .2, .5. 1.1. 1.3. 2.0. .3. 0.0, -.1.0,0, .3 15 6 2 - 
.1,0,0, .1. .4. .3,9 
96 DATA 0. .2. -.1.0. - 

-1, .1, .l,.l. .1, .2,0,0.0.0,0,0,0,0,0.0.0.0.0.0.0.0,0.0 
98 DATA 0,0.0,0.0,0.0.0,0.0.0,0.0,0.0.0,0 0 0 0- 
.1.0,0, .7, .7. .6. .1,0. .1 

100 DATA 0. .1, .2. .3, .6. .8. 1 .0, .2.0.0. - 
.1.0.0. . 1. .1$, .3. .1.0.0.0.0. .2. .2. .5 
102 DATA 0, .1, - 

.1.0.0.0. .1. .1, .1, . 1.0.0.0.0,0,0.0.0.0.0.0.0.0.0,0,0.0,0.0,0 
104 DATA 

0.0,0,0,0,0,0,0,0.0,0,0,0.0.0,0,0,0.0,0.0, .5, .4, .4. .1,0,0.0 0 1 1 
106 DATA .4, .5, .6. .1,0.0,0.0,0,0. .05. .1. 0.0.0.0.0, . 1 /. l' 301 - 

,1.0,0,0.0,0 

108 DATA -1,0,0,0,0,0.0.0,0.0.0.0,0,0,0,0,0 0,0 0 0 0 

120 DIM V(9.9,9) 

130 FOR Z-X TO 9 STEP 1 

140 FOR Q-l TO 9 STEP 1 

145 FOR 'J=l TO 9 STEP I 

147 READ V(Z,Q,U) 

222 NEXT U 

223 NEXT 0 

224 NEXT Z 

225 PRINT II. * VOLTAGE DATA* 

226 DIM VO<9,9> 

227 FOR Q=l TO 9 STEP 1 

228 FOR R=l TO 9 STEP 1 
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229 FOR U-l TO 9 STEP 1 

230 LET H=V<U,R,Q) 

231 LET VO(U,R}=H 

232 NEXT U 

233 NEXT R 

234 MAT PRINT #1, VO , 

235 NEXT 0 

254 ! FFT THE ROWS OF SAMPLED VOLTACES IN THE X DIRECTION 

255 DIM KR(9) 

256 DIM KI(9) 

257 DIM R(9) 

258 DIM RV(9,9,9> 

259 DIM IV(9.9.9) 

260 FOR V=0 TO 9 STEP 1 
270 FOR M*l TO 9 STEP 1 
280 FOR N=l TO 9 STEP 1 
290 LET Y*V<M.N,V) 

300 LET R(N)=Y 
305 NEXT N 

310 EXTERNAL SUa FFT ( DIM { ) , DIM {) , D JM{ ) ) 

320 CALL FFT IR() , MR ( ) , M I ( ) ) 

330 FOR N-l TO 9 STEP 1 

340 LET Y-MR (N} 

350 LET RVIM.N ( V|=Y 

360 NEXT N 

370 FOR M=l TO 9 STEP I 

380 LET G=MI(N) 

390 LET IV(M,N.V) 

400 NEXT N 

410 NEXT K 

415 NEXT V 

417 DIM RVO(9,9) 

419 DIM IVO(9.9) 

420 PRINT #1. *RV* 

421 FOR 0*1 TO 9 STEP 1 
42 2 FOR R=l TO 9 STEP 1 

423 FOR U=l TO 9 STEP 1 

424 LET IURV(U.R,Q) 

425 LET RVO(U.R)=H 

426 NEXT y 

427 NEXT R 

428 KAT PRINT »1. RVO, 

429 NEXT 0 

430 PRINT fll, "IV- 

431 FOR 0=1 TO 9 STEP 1 

432 FOR R=l TO 9 STEP 1 

433 FOR U=l TO 9 STEP 1 

434 LET H=IV(U,R,Q) 

435 LET IVO(U,RUH 

436 NEXT U 
4 37 NEXT R 

438 MAT PRINT HI. IVO. 

439 NEXT 0 

tit Dr" R Rvn C ?T S ° K ™ E SAMPLE ° WLTAGE IN ™ E 2 DIRECTIOM 

470 DIM LRVC9.9.9) 

475 FOR Vil TO 9 STEP 1 

480 FOR N=l TO 9 STEP I 
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490 FOR H = l TO 9 STEP L 
500 LET Y=RV(M.N,V) 
S10 LET R(M)=Y 
520 NEXT M 

530 CALL FFT(R(),MR(),MIO) 

540 FOR M= 1 TO 9 STEP 1 

550 LET H-MR (M) 

560 LET RRV(M,N,V)=H 

570 NEXT M 

580 FOR M-l TO 9 STEP I 

590 LET G*MI (H) 

600 LET IRV<M.N,V)=G 

610 NEXT M 

620 NEXT N 

625 NEXT V 

626 DIM RRVO(9,9) 

627 DIM IRV0<9, 9) 

630 PRINT II, * RRV" 

631 FOR 0=1 TO 9 STEP 1 
63 2 FOR R=l TO 9 STEP I 
63 3 FOR U =1 TO 9 STEP 1 
634 LET H=RRV(U.R.QJ 

637 LET RRVO<U,R)=H 

638 NEXT U 
63 9 NEXT R 

640 MAT PRINT SI, RRVO, 

641 NEXT Q 

650 PRINT HI . ' IRV- 

651 FOR O^l TO 9 STEP I 

653 FOR R=?l TO 9 STEP 1 

654 FOR U=l TO 9 STEP I 

655 LET H- IRV(U,R,OI 

656 LET lRVO{U,R)=H 

657 NEXT V 

658 NEXT R 

659 MAT PRINT HI. IRVO. 

660 NEXT O 

661 DIM RIV(9,9.9> 

662 DIM I1V(9,9,9> 

685 FOR V=l TO 9 STEP 1 
690 FOR N=l TO 9 STEP I 
700 FOR H=l TO 9 STEP 1 
710 LET Y=IV|M,N,V) 
760 LET R(M)=Y 
770 NEXT M 

830 CALL FFT { R () , MR ( ) , MI { ) j 

840 FOR M=l TO 9 STEP 1 

850 LET H=MR (M) 

860 LET RIV(M,N,V}*H 

870 NEXT M 

872 FOR M= 1 TO 9 STEP I 

873 LET G=MI(M> 

874 LET IIV(M,N,V)=G 

875 NEXT M 

876 NEXT N 

877 NEXT V 

878 DIM RIVOC9. 9) 

879 DIM IIVO<9.9) 
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880 PRINT II, -Riv 

881 FOR 0=1 TO 9 STEP 1 

882 FOR R=l TO 9 STEP 1 

883 FOR U=l TO 9 STEP I 

884 LET H=RIV(U,R,Q) 

885 LET RlVO(a,RJ=H 

886 NEXT U 

887 NEXT R 

888 MAT PRINT II. RIVO, 

889 NEXT Q 

890 PRINT II. -IIV 

891 FOR Q=\ TO 9 STEP 1 

892 FOR R=l TO 9 STEP 1 

893 FOR U:l TO 9 STEP 1 

894 LET H=IIV(U,R,Q) 

895 LET IIVO(V.R)=H 

896 NEXT 1/ 

897 NEXT R 

898 MAT PRINT II. IIVO. 

899 NEXT Q 

900 DIM RVA(9,9,9| 

901 DIM IVA(9,9,9J 
904 DIM RVAO<9. 9) 
906 DIM IVAO(9.9) 

908 FOR Q= 1 TO 9 STEP 1 

910 FOR R^l TO 9 STEP 1 

911 FOR \J~l TO 9 STEP 1 

912 LET H-IIV<U,R,Q> 

913 HM-1) *H 

914 LET G=RRV(U,R,Q> 

915 LET L=G+H 

916 LET RVA(U,R,Q)=L 

917 NEXT U 

918 NEXT R 

919 NEXT Q 

920 FOR 0-1 TO 9 STEP 1 

921 FOR R = l TO 9 STEP 1 

922 FOR U - I TO 9 STEP 1 

923 LET H=IRV(U,R,QJ 

924 LET U-RXV(U.R.O) 

925 LET IVA(U,R,Q}=H*I. 

927 NEXT U 

928 NEXT R 
930 NEXT 0 

932 PRINT #1. -RVA* 
934 FOR Q-i TO 9 STEP 1 
936 FOR R=l TO 9 STEP 1 
938 FOR U=l TO 9 STEP 1 
940 LET H=RVA<U,R,Q) 
942 LET RVAOCU.R)=H 
944 NEXT V 
946 NEXT R 

948 MAT PRINT #1, RVAO. 
950 NEXT Q 
952 PRINT 11, -IVA- 
954 FOR 0=1 TO 9 STEP 1 
956 FOR R^l TO 9 STEP 1 
958 FOR 0=1 TO 9 STEP i 
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960 LET H=IVA(U,R,Q} 
962 LET 1VA0(U.R)=H 
964 NEXT U 
966 NEXT R 

968 MAT PRINT II, XVAO, 
970 NEXT Q 

llVoil RVAyS W 9,9? SA " ,PLED V0LTAGES IN ™ E Y DIRECTION 

990 DIM IRVAY<9,9.9) 
992 DIM RVAYO(9,9) 
994 DIM IRVAYO{9,9) 

1012 FOR MM TO 9 STEP 3 

1013 FOR N=l TO 9 STEP I 

1014 FOR V*l TO 9 STEP I 

1015 LET Y=RVA(M.N,V) 

1016 LET R(V)=Y 

1017 NEXT V 

1018 CALL FFT<R<},MRU.KIO) 

1019 FOR V=l TO 9 STEP 1 

1020 LET H=HR(V> 

1021 LET RVAY (H,N,V>eH 

1022 NEXT V 

102 3 FOR V=l TO 9 STEP 1 

1024 LET G-MI(V) 

1025 LET IRVAY(M. N, V) *Q 

1026 NEXT V 

1027 NEXT N 

1028 NEXT M 

1030 PRINT 81, -RVAY" 
1032 FOR 0=1 TO 9 STEP 1 
1034 FOR R=l TO 9 STEP 1 
1036 FOR U=l TO 9 STEP 1 
1038 LET H=RVAY(U,R.0) 
1040 LET RVAYOtU.RJ^H 
1042 NEXT U 
1044 NEXT R 

1046 MAT PRINT *l,RVAYO. 
1048 NEXT Q 
1050 PRINT HI. * I RVAY" 
1052 FOR 0-1 TO 9 STEP 1 
1054 FOR R=l TO 9 STEP 1 
10S6 FOR U=l TO 9 STEP 1 
1058 LET H=IRVAY(U,R,Q) 
1060 LET IRVAYO(U.R)=H 
1062 NEXT U 
1064 NEXT R 

1066 MAT PRINT *l,IRVAYO, 
1068 NEXT 0 
1070 DIM RIVY<9,9,9) 
1080 DIM IIVY{9.9. 9) 

1085 DIM RIVYO<9,9) 

1086 DIM IIVYO<9,9) 
1090 FOR M=l TO 9 STEP 1 
1100 FOR N=l TO 9 STEP 1 

1138 FOR V=l TO 9 STEP 1 

1139 LET Y=IVA(M.N,V> 

1140 LET R(V)=Y 

1141 NEXT V 
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1142 CALL FFTIRO.HRO.Mim 

1143 FOR V=l TO 9 STEP I 

1144 LET H*KR(V) 

1145 LET RIVY{M,N,V)=H 

1146 NEXT V 

1147 FOR V=l TO 9 STEP 1 

1148 LET G=MI(V) 

1149 LET HVY(M,N,V)=G 

1150 NEXT V 

1151 NEXT N 

1152 NEXT M 

1153 PRINT HI, -RIVY" 
1160 FOR 0=1 TO 9 STEP 1 
1162 FOR R=l TO 9 STEP 1 
1164 FOR U=l TO 9 STEP 1 
1166 LET H«RIVY(U,R,Q) 
1170 LET RIVYO(U.R)=H 
1172 NEXT U 

1174 NEXT R 

1178 MAT PRINT II, RIVYO. 
U80 NEXT 0 
1185 PRINT 8i . • IIVY* 
1190 FOR Q=l TO 9 STEP 1 
1200 FOR R=l TO $ STEP 1 
1210 FOR U=l TO 9 STEP 1 
1212 LET H= IIVY<U,R.Q) 
1214 LET irVYO(U.R>=H 
121$ NEXT U 
1218 NEXT R 

122 0 HAT PRINT »l,IIVYO, 

12 2 2 NEXT 0 

1230 DIM YRVA(9,9,9) 

1240 DIM YIVM9.9, 9) 

1241 FOR 0=1 TO 9 STEP 1 
1250 FOR R=l TO 9 STEP I 
1260 FOR U=l TO 9 STEP i 
1270 LET L=IIVY<U,R,Q> 
1280 LET B=RVAY(U.R.Q) 
1290 LET YPVA(U.R.O)=B-L 
1300 NEXT U 

1310 NEXT R 
1320 NEXT 0 

1330 FOR 0^1 TO 9 STEP 1 
1340 FOR R=I TO 9 STEP 1 
134S FOR Url TO 9 STEP 1 
1350 LET H=RIVY<U,R,Q) 
1360 LET L=IRVAY(U.R.Q) 
1370 LET YtVA<U.R,Q> =L*H 
1380 NEXT U 
1390 NEXT R 
1400 NEXT 0 
1410 PRINT ¥1, ' YRVA " 
1412 DIM YRVAO(9,9) 
1414 DIM YIVAO(9,9) 
1420 FOR 0=1 TO 9 STEP 1 
1430 FOR R*l TO 9 STEP 1 
14 40 FOR U=i TO 9 STEP 1 
1450 LET H=YRVA(U.R.Q> 
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1460 LET YRVAO<U.R)=H 
1470 NEXT U 
W80 NEXT R 

1490 MAT PRINT »l,YRVAO, 
1500 NEXT Q 
1510 PRINT #1. "YIVA- 
1520 FOR Q=l TO 9 STEP 1 
1530 FOR M TO 9 STEP 1 
1540 FOR 0=1 TO 9 STEP 1 
1545 LET H=YIVA(U,R,0) 
1550 LET YIVAO<U,R>=H 
1560 NEXT U 
1570 NEXT R 

1580 MAT PRINT *l,YIVAO, 
1590 NEXT Q 

Ittl 122^12 ™ E DISCRETE SPECTRUM OF THE SYSTEM FUNCTION AND THE SXNC 

llltJr FUNCTI0N 0F THE SAMPLE HOLD CORRESPONDING TO THE FINITE 

DLTECTOR 

2057 ! DIMENSIONS 

2060 DIM SFH(9,9,9) 

2061 DIM SFHO{9.9) 

2062 DIM SINC (9,9.9) 

2063 DIM SINCO(9.9) 

2065 LET SFH{4. 4,4)=4*PI 

2066 FOR M=-4 TO 4 STEP 1 
2070 FOR V=-4 TO 4 STEP 1 

2080 FOR N^-4 TO 4 STEP 1 

2081 H=.04 

2082 J=ABS(MJ +ABS1N) «-ABS(V) 
2085 IF J=0 THEN CO TO 2098 

2090 T*4*PZM(2*PI-N/9*1/.1)~2«((2*PI«V/9M/ 1>~2)) 

2092 B-(2'PI*N/9M/.1)-2M2*PI*M/9*1/C)-2*(2*PI-V/9M/ 1>~2 

2095 LET SFH<M*5,N*S, V+SJ^T/B 

2098 C=ABS<N)*ABS{V) 

2100 IF G^0 THEN GO TO 2114 

2101 A=l 

2105 IF N=0 THEN GO TO 2107 

2106 A=SIN<2*PI*10*N/9» . 1 ) / <PI *N/9 * 10) 

2107 3-1 

2108 IF V=0 THEN CO TO 2110 

2109 B=SIN(2*PIU0*V/9*.U/ (Pr*V/9*l0> 

2110 H=A*B 

2111 IF N=0 THEN H=.2'H 

2112 IF V=0 THEN H=,2*H 
2114 LET SINC<M*S,N«S,V*5>=H 

2130 NEXT N 

2131 NEXT V 

2132 NEXT M 

2135 PRINT tfl, ■ SFH" 

2136 FOR 0=1 TO 9 STEP 1 

2137 FOR R-l TO 9 STEP 1 

2138 FOR 11= 1 TO 9 STEP 1 

2139 LET S=SFH(U < R.O) 

2140 LET SFHO(U.R)*S 

2141 NEXT U 

2142 NEXT R 

2143 MAT PRINT «1. SFHO, 

2144 NEXT Q 
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2145 PRINT *1, 

2146 FOR 0=1 TO 9 STEP I 

2147 FOR R=l TO 9 STEP 1 

2148 FOR U^l TO 9 STEP 1 

2149 LET S~SINC(U,R,Q) 

2150 LET SINCO{U,R)=S 

2151 NEXT U 

2152 NEXT R 

2153 HAT PRINT §1, SINCO, 

2154 NEXT 0 

llti Sno7 ° ISCRETE SPECTRUM 0F THE SYSTEM FUNCTI °" ™> ™ sinc 

2160 DIM HR<9.9,9) 

2170 DIM HRO(9.9) 

2171 DIM SrNCR(9.9,9> 
217S FOR V=l TO 9 STEP 1 
2180 FOR M-l TO 9 STEP 1 
2190 FOR N=l TO 9 STEP 1 

2200 LET Y^SFH(M,N,V) 

2201 LET H=SiNC(M,N,V> 

2210 IF Y-0 THEN GO TO 2212 

2211 y=l/y 

2212 IF H=0 THEN GO TO 2250 
2221 H=l/H 

2230 LET HR(M,N.V)=Y 
223S LET SINCR<K.N,V)=H 
2240 NEXT N 
2250 NEXT M 
2260 NEXT V 

2310 PRINT »1, " HR * 

2311 FOR 0^1 TO 9 STEP 1 

2312 FOR R=l TO 9 STEP 1 

2313 FOR U= 1 TO 9 STEP 1 

2314 LET S*HR<U,R,Q) 

2315 LET HRO(U>R3=S 

2316 NEXT U 

2317 NEXT R 

2320 MAT PRINT #l,HRO. 

2321 NEXT Q 

FTOCTION VIDE ^ TRANSP0RJ<ED DATA BY THE TRANSFORM OF THZ SYSTEM 

2345 ! AND THE SINC FUNCTION 

3030 DIM FR<9.9,9) 

3050 DIM FM9.9.9} 

3052 DIM FRO<9,9) 

3054 DIM FIO(9,9> 

3056 FOR V=l TO 9 STEP I 

3057 FOR M=l TO 9 STEP 1 

3058 FOR N=i TO 9 STEP 1 

3059 T=YRVA (M, N. V) 

3060 S^HR(M.N.V) 
3062 L^SINCR(K.N,VJ 

3065 K=S*T*L 

3066 LET FR<M,N,V)=K 

3067 NEXT N 

3068 NEXT M 
3070 NEXT V 

3080 FOR Vri TO 9 STEP 1 
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3081 FOR M=i TO 9 STEP L 

3082 FOR N=l TO 9 STEP i 

3083 H=SINCR(M,N. V} 

3093 K= HR(M,N,V) 

3094 L= YIVA<M,N,V) 
309$ £=K*L*H 

3096 LET FICM r N.V)=E 

3097 NEXT N 

3098 NEXT M 
3100 NEXT V 

3130 PRINT II. -FR* 

3131 FOR 0=1 TO 9 STEP 1 

3132 FOR R=l TO 9 STEP 1 

3133 FOR V-l TO 9 STEP i 

3134 LET 11 = FR(U.R,Q> 

3135 LET FRO(U.R>=H 

3136 NEXT U 

3137 NEXT R 

3138 MAT PRINT #1, FRO, 

3140 NEXT Q 

3141 PRINT fl, -FI" 

3142 FOR 0=1 TO 9 STEP 1 
314 3 FOR R = 1 TO 9 STEP 2 

3144 FOR U=i TO 9 STEP 1 

3145 LET H=FI(U,R r Q) 

3146 LET FIO<U.R)*H 
3 156 NEXT U 

3157 NEXT R 

3158 MAT PRINT fl. FIO, 
3160 NEXT Q 

3170 ! INVERSE TRANSFORM THE ROWS IN THE X DIRECTION 
3180 DIM RF(9.9,9> 
3185 DIM RFO{9,9) 
3187 DIM MlFO(9,9) 
3190 DIM MIF(9,9.9] 

3195 DIM 1(9) 

3196 FOR V=l TO 9 STEP 1 
3 200 FOR M= 1 TO 9 STEP 1 
3210 FOR N=l TO 9 STEP 1 
3220 LET Y=FR(M,N,V> 
3230 LET R(N)^Y 

3 240 NEXT H 

32S0 FOR N=l TO 9 STEP 1 
3260 LET Y=FI<M,N,V> 
3270 LET I(NMY 
3280 NEXT N 

328S EXTERNAL SUB IFT [ DIM ( ) , DIM { ) , DIM ( ) OIM{)) 
3290 CALL IFT(R ( ) , I { ) , MR ( ) . MI ( ) ) 
3300 FOR N=l TO 9 STEP 1 
3310 LET Y=KR (N) 
3320 LET RF(M.N,V}=Y 

333 0 NEXT N 

334 0 FOR N=l TO 9 STEP 1 
33S0 LET Y^MKN) 

3360 LET HIF(M,N.V)=Y 
3 370 NEXT N 
3380 NEXT M 
3385 NEXT V 
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3390 PRINT tl, -RF- 

3391 FOR Q=i TO 9 STEP I 

3392 FOR R=a TO 9 STEP 1 

3393 FOR U=l TO 9 STEP I 
3 3394 LET H=RF(U,R,Q> 

3395 LET RFO(U.R)=H 

3396 NEXT U 

3397 NEXT R 

3400 MAT PRINT »1, RFO, 
*U 3405 NEXT Q 

3410 PRINT II,- MIF • 

3411 FOR Q=l TO 9 STEP 1 

3412 FOR R=l TO 9 STEP 1 

3413 FOR U=l TO 9 STEP 1 

3414 LET H. MIF(CI,R.Q> 

3415 LET MIFO{U,R^H 

3416 NEXT 0 

3417 NEXT R 

3420 MAT PRINT 11, M IFO 
-'0 3421 NEXT O 

3422IINVERSE TRANSFORM THE ROWS IN THE V DI R.ECTTf»M 
3430 DIM RFY(9.9,9) DIRECTION 

3432 DIM MIFY(9,9.9) 

3433 DIM RFYOf9,9> 
-5 343 4 DIM MIFYO<9,9} 

3435 FOR M=l TO 9 STEP I 
3440 FOR N=a TO 9 STEP 1 
3450 FOR V=l TO 9 STEP 1 
3460 LET Y=RF <M, N, V) 
->t> 3470 LET R(V)^Y 
3480 NEXT V 

3490 FOR V=l TO 9 STEP 1 
3500 LET Y=MIF<M,N.V) 
3510 LET ICV)=Y 
3520 NEXT V 

3525 EXTERNAL SUB I FT (DIM O . DIM { ) DIMM DIMM* 

3526 CALL IFT<R ( > . r U , MR U , MI < >) "° > 

3527 FOR V-l TO 9 STEP 1 

3528 LET Y*MR<V> 

40 3529 LET RFY(M.N,V)= Y 

3530 NEXT V 

3531 FOR V= I TO 9 STEP 1 

3532 LET Y-MI(V) 

3533 LET MIFY(M,N.V)=Y 

3534 NEXT V 
353 5 NEXT N 

3536 NEXT M 

3537 PRINT #1. -RFY- 

35 38 FOR Q=l TO 9 STEP 1 
->U 3539 FOR R=l TO 9 -STEP 1 

3540 FOR V-~l TO 9 STEP 1 

3541 LET H=RFY(U,R,Q) 

3542 LET RFYO(U,R)=H 

3543 NEXT U 
55 3S44 NEXT R 

3545 MAT PRINT 81. RFYO 

3546 NEXT 0 

3547 PRINT ftl, -MIFY* 
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3550 FOR Q=l TO 9 STEP 1 
355S FOR R*l TO 9 STEP 1 
3560 FOR U=l TO 9 STEP 1 

3565 LET H^MIFYCU, R, Q) 

3566 LET MIFYO(U,R)=H 

3567 NEXT U 

3568 NEXT R 

3570 MAT PRINT !l,MIFYO, 
3575 NEXT Q 

1^ TT SF0RM ™ E C0LUMNS IN ™ E 2 DIRECTION 

3590 DIM FO(9,9> 

3592 FOR VxJ TO 9 STEP 1 

3593 FOR N-l TO 9 STEP 1 

3594 FOR M=l TO 9 STEP 1 
3600 LET Y=RFY(M,N,V) 

4546 LET R(M>=Y 

4547 NEXT M 

4548 FOR M=l TO 9 STEP 1 

4549 LET Y=HIFY Hi. N, V) 

4550 LET I(K)=Y 

4552 NEXT K 

4553 EXTERNAL SUB IFTZ (DIM ( ) . DIM f ) , DIM ( ) DIMO) 

4554 CALL IFTZ (R < > , I O , MR ( ) . MI {) ) 

4556 FOR K-l TO 9 STEP I 

4557 LET Y=MR (M) 
4559 LET F(M,N,V)=Y 
4600 NEXT M 

4602 NEXT N 

4604 NEXT V 

4605 PRINT *1, "F* 

4610 FOR 0=1 TO 9 STEP 1 
4620 FOR R=l TO 9 STEP 1 
4630 FOR U=l TO 9 STEP I 
4635 LET H=F(U.R,Q) 
4640 LET FO(U, R)=M 
4650 NEXT U 
4660 NEXT R 
4666 MAT PRINT 81,FO. 
4670 NEXT 0 

till DIH^.9?;, THE U(Z1 CONVOL0TI0N ^ NORMALIZE THE RECONSTRUCTION 
4780 DIM CFX3(9.9> 

4809 FOR V=l TO 9 STEP 1 

4810 FOR N=l TO 9 STEP t 
4871 LET K=F(1.N,V> 
4873 K=K/10E8 

4974 LET CF(1,N,V)=K 
5915 NEXT N 

6616 NEXT V 

6617 FOR V=l TO 9 STEP 1 
6620 FOR M=2 TO 9 STEP 1 
6630 FOR N=l TO 9 STEP 1 
6640 LET K=F(M,N.V) 

6642 P=0.5*PI/<M-l)^2/C-2 
6650 K=K/P 
6660 LET CF(M,N.V)=K 
6670 NEXT N 
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6680 NEXT M 

6681 NEXT V 

6690 PRINT »X. " RECONSTRUCT ION * 

6691 FOR Q^l TO 9 STEP 1 

6692 FOR R=l TO 9 STEP I 

6693 FOR U=l TO 9 STEP 1 

6694 LET H=CF(U,R,Q) 

6695 H=H/CF<5,5,6) 

6696 LET CFO(U,R)=H 

6697 NEXT V 

6698 NEXT R 

6700 MAT PRINT ll.CFO. 
6705 NEXT Q 
6710 END 

6722 SUB FFT ( R { ) . MR O . MI ( ) ) 
6726 FOR M=l TO 9 STEP 1 
6730 A=0 

6740 FOR N=l TO 9 STEP 1 
6750 LET H=R(N) 

5760 B-K*COs;2*PIMM-S) MN-5)/9) 

6770 A-A+B 

6780 NEXT N 

6790 A=A/9 

6800 LET MR<M>*A 

5810 NEXT M 

6820 FOR M=l TO 9 STEP t 
6830 A-0 

6840 FOR N = l TO 9 STEP 1 
6880 LET H=R(N> 
6890 H=-H 

6900 B=H*SIN(2*PIMM-5)MN-5)/9) 

6910 A=A+B 

6920 NEXT N 

6930 A-A/9 

6935 LET MKMJ-A 

6940 NEXT M 

6950 END SOB 

6960 SUB XFT(RU.XO # MRO.MX<)) 

5970 DIM MRRI9) 

6980 DIM MRI19) 

6990 DIM MIR{9) 

7000 DIM Mil (9! 

7010 FOR N=l TO 9 STEP 1 

7020 A=0 

7030 FOR M*l TO 9 STEP 1 
7040 LET G-R (M) 

7050 B=G*COS(2*PIMM-5)MN-5)/9> 

7060 A=A+B 

7070 NEXT M 

7080 LET MRR(N) =A 

7090 NEXT N 

7100 FOR N=l TO 9 STEP 1 
7110 A=0 

7120 FOR M=l TO 9 STEP 1 
7130 LET G=R(M| 

7140 0=G»SIN(2»PIMM-5) MN-5) /9) 
7150 A^AtB 
7160 NEXT K 
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7170 LET MXI(N)-A 
7180 NEXT N 

7190 FOR N=i TO 9 STEP 1 
7200 A=0 

7210 FOR M=l TO 9 STEP 1 
7220 LET G=I(M) 

7230 B=G*COS(2*PIMM-S) MN-S)/9> 

7240 A-A+B 

7250 NEXT M 

7260 LET MIR<N)=A 

7270 NEXT N 

72B0 FOR N=l TO 9 STEP 1 
7290 A^0 

7300 FOR M=L TO 9 STEP 1 
7310 LET G^I(M) 

7320 B*C*SlN(2*PIMM-5)*(N-5)/9) 

7330 A=A*8 

7340 NEXT M 

7350 LET HII<N)=A 

7360 NEXT N 

7365 MAT MIX=(-ll-MIJ 

7375 MAT MR=MRR+MI I 

7385 MAT MI=MIR*MRI 

7400 END SUB 

7410 SUB IFTMRO.IO.MRO.Mim 

7420 DIM MRR ( 9 ) 

7430 DIM MRI(9> 

74 40 DIM MIR (9) 

7450 DIM MU<9> 

7460 FOR N=l TO 9 STEP 1 

7470 A=0 

74 80 FOR M-l TO 9 STEP 1 
7490 LET G-R (M) 

7500 B=G*COS(2*PIMM-5| *(-N* l)/9) 

7510 A=A+B 

7520 NEXT M 

7S30 LET MRR (N) =A 

754 0 NEXT N 

7550 FOR N=l TO 9 STEP J 
7^60 A=0 

7570 FOR M=l TO 9 STEP 1 
7580 LET G=R<M) 

7590 B=G*SIN(2*PIMM-5) • (-N-D/9) 

7600 A«=A+B 

7610 NEXT M 

7620 LET MRI(NJ=A 

7630 NEXT K 

7640 FOR N=l TO 9 STEP 1 

7641 A=0 

7642 FOR M*l TO 9 STEP 1 
7650 LET C=I<M) 

7660 B=C*COS(2»PI* (M-S) • (-N+1) f9) 

7670 A=A+B 

7680 NEXT M 

7690 LET MIR<N)=A 

7700 NEXT N 

7710 FOR N=l TO 9 STEP 1 
7720 A=0 
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77J0 FOR M=l TO 9 STEP 1 
7740 LET G=I (M) 

7750 B=C*SINt2*PIMM-SI*(-N*l»/9) 
7760 A=A*B 
5 7770 NEXT M 

7780 LET MII(H)=A 
7790 NEXT N 
779S MAT MII={-1) »MH 
7800 MAT MR=MRR*MI I 
10 7810 MAT MI=MRItMIR 
7880 END SUB 
5 



In an embodiment, the matrix inversion and/or the reiterative 

aSSSH a r 0 used ,n , combina,ion whh the Fou '- Trl£T 

se!Z£Jl<! : h Xam p ple - m < hc casc of NM* measurement of the 

a serfs of fi « ^ TranSf ° rm AIg0ri,hm is uscd «° «^ula«e 

20 detected by usin \ . ^« 

»s.gncd to a voxel or the value of the derivative of the magn L 
suscepuMiiy as a function of voxel position, respectively tE 
roagnet.c susceptibility of the voxcl( S ) corresponding to the peaks or 

25 L^Tre^nT^s <° «* 

^ 4UC,,C ' UU-^J)) and the remaining voxels are 

elr is " ?" $ ' gnal ^ C ° mparCd *° ,he aC,ual si ^» die 

error .s determined. The error is back projected over the voxels using 

30 Reite^vTZ"" 5 ^ Ma,rix fnvcrsio " Algorithm or 8 

3encv Th IT*" 6 " 0 susce P ,ibili <y for each Larmor 

frequency. The procedure is rciteratively until it stabilizes This 

35 1 > d ,nde P endenl 'y <via by the Larmor frequency) from the 
35 pat.al reconstructs (via the geometric system function . the 
recon Slru tto d , versus ^ ^ tcrms J 

,„ .K ,' M h , e SU P er,orit y of digital over analogue is generally known 
to those skilled in the art of signal processing. 

40 The Nyquist Theorem With The Determination Of The Spatial 

Resolution 

The derivation of Eq. (30) demonstrates that the system 
45 suscept,b.l,ty funcon. It is well known in the art of signal processing 
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that if such a filter passes all frequencies for which an input function 
has significant energy, then ihe input function can be recovered 
completely from samples of the filtered function taken at the Nyquist 
rate. This premise embodies the Nyquist Sampling Theorem The 
spectrum of the system function (Eq. (2)) is derived in APPENDIX V 
and shown in FIGURE lc. This function is a band-pass for all 
frequencies of the magnetic susceptibility function where k r and *. 
are comparable. Thus, the magnetic susceptibility function can be 
recovered by sampling the continuous voltage function given by Eq 
(29) at the Nyquist rate, twice the highest frequency of the magnetic 
susceptibility function, in each spatial dimension over the sample 
space for which the function has appreciable energy. Sampling 
operations other than the present operation and the negligible error 
encountered by not sampling over the entire sample space are 
15 discussed in McC. Siebert [II] and the references therein disclosed 
which are all incorporated herein by reference. In the absence of 
noise, the spectrum of the magnetic susceptibility function can be 
completely recovered if the detector spacing frequency is equal to the 
Nyquist rate which is twice Ihe highest frequency of the magnetic 
20 susceptibility function, and this represents the limit of resolution 

However, the density of ihe detector spacing is limited by noise The 
three-dimensional magnetic susceptibility map is a reconstruction 
from independent recordings at independent detector spatial locations 
relative to the voxels of the image space where two detector signals 
are independent if they are sufficiently spatially displaced from each 
other such that the difference in signal between the two detectors due 
to a test voxel is above the noise level. The resolution based on 
signal-to-notse arguments is discussed in the Contrast and Limits of 
Resolution Section. 



25 



30 



40 



Contrast and Limits of Resolution 



The ability to visualize a structure in a noise-free 
environment depends, among other factors, on the local contrast C 
35 which is defined as 



~ / (63) 



where / is the average background intensity and A/ is the intensity 
variation in the region of interest. The contrast for ReMSI is typically 
greater than 20%, which compares favorably with computed 
tomography where the conirast of soft tissue is about 0.1%. The main 
source of NMR image (also called magnetic resonance images (MRI)) 
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contrast is T x and 7\ which depend on tissue types. In an 
embodiment. T, and 7\ pulse sequences are provided by magnetic 
field gradient coils and are used to alter the intensity of the 
secondary RF field as a function of time to provide image contrast. 
5 The contrast may be increased by using RF pulse sequences to 

polarize the protons which encode T % (spin/lattice) and 7, (spin/spin) 
relaxation time information in the data of the secondary RF field in a 
manner straightforward to those skilled in the NMR art to apply the 
fields and detectors described for this invention. In an embodiment, 
10 the intensity of the bulk magnetization M of each voxel having the 

corresponding magnetic moment is plotted for each Larmor frequency 
rather than the magnetic susceptibility wherein T x (spin/lattice) and 
7\ (spin/spin) relaxation times encode contrast via intensity 
variations that arc tissue specific and dependent on r, and/or 7\. 
15 Contrast, however, is not a fundamental limit on visualization 

since it can be artificially enhanced by. for example, subtracting part 
of the background or raising the intensity pattern to some power. 
Noise, however, represents a fundamental limitation on the abiiitv to 
visualize structures. The signai~io-noise ratio, a basic measure of 
visualization, is determined by the ratio of the desired intensity 
variations to the random intensity variations, which are governed by 
the statistical properties of the system. The signal-to-noise ratio 
(SNR) is defined as 

SNR-f-cl (64) 

25 where 0, is the standard deviation of the background intensity 

representing the rms value of the intensity fluctuations. The noise 
properties of the RcMSI phenomenon involve additive noise only 
principally from thermal (i) noise in the RF measurement circuits of 

the antennas of the detector array and from the fluctuations of the 

30 primary magnetic field. 

It can be seen from Eq. (23) that the total magnetic 
susceptibility map can be determined from the total magnetic 
moment map if the magnetizing flux map is known. Since 
**l<r , -IO" s «l, the magnetic susceptibility x causes a small 

35 perturbation to the magnetizing flux. Thus, the magnetizing flux must 
be stable to at least the magnitude of the signal contribution due to 
the magnetic susceptibility. A feature of superconducting magnets is 
their extreme stability. The object to be imaged is magnetized with a 
highly stable magnet such as a superconducting magnet? In this case. 

40 a magnetic field stability of 10"* % over a month's time is feasible 
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The magnetizing means can also possess a means to add a component 
of modulation to the magnetizing field at frequencies below those 
which would induce eddy currents in the tissue which would 
contribute significant noise to the secondary magnetic field signal. 
5 Such modulation would cause an in-phasc modulation of the RF 

secondary magnetic field signal which can be detected and filtered to 
improve the signal to noise ratio, or the modulation over space can be 
used to encode spatial information. 

It can further be seen from Eqs. (20-23) that the total magnetic 

10 moment map can be determined from the NMR frequencies as input. 
Since x=* W~* - W s « I- the magnetic susceptibility x causes a small 
perturbation to each NMR frequency (Eqs. (20-23)). However, 
extremely high frequency resolution is routine. The frequency 
resolution of RF signals exceeds 10"* %. Small antennas measure the 

15 RF signals as point samples without significant decrease in the signal 
to noise ratio relative to large antennas by using impedance matching 
while minimizing resistive losses by using superconducting reactance 
elements, for example. In an embodiment, cross talk between 
antennas is ameliorated or eliminated by time multiplexing the signal 

20 detection over the array of antennas. External sources of RF noise can 
be ameliorated by placing the RcMSI scanner in a shielded room 
(Faraday cage). The quality of the image (i.e. the signal to noise ratio 
of the image) can be increased by repeating the reconstruction over 
multiple time points wherein each data set of a given time point is the 

25 set of matrices of the intensity variation over the sample space of the 
RF field of the bulk magnetization M of each voxel having the 
corresponding magnetic moment at that synchronized (common) time 
point. Sources of detector error are random, and the noise averages 
out as the number of detectors increases. Typically, the noise is 

30 suppressed by a factor of the inverse square root of the number of 
detectors. Thus, the effective limit of fiux detection is increased by a 
factor of the square root of the number of detectors for a constant 
SNR. 

The resolution depends on the extent that the field of the ring, 
35 shell, or sphere of dipoles differs from that of a single dipole at the * 
center. The plot of the three cases of the field of a ring, shelK and a 
sphere or dipoles each of radius R and magnetic moment m = I0 4 Gem- 
given by Eq. (1.14). Eq. (11.17), and Eq. (IV.16) of APPENDIX I 
APPENDIX II, and APPENDIX IV, respectively, as a function of radius 
40 R where the position of the center of the ring, shell, and sphere 

relative to the detector is the point (0.0.10) is given in FIGURES 2, 4. 
and 6, respectively. The plot of the three cases of the field of a ring, 
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20 



she! and sphere of dipoles of radius R = 0.2cm and magnetic moment 
m = IO Gem given by Eq. (1. 14), Eq. (11.17), and Eq. (IV 16) of 
APPENDIX i. APPENDIX H, and APPENDIX IV, respectively as a 
function of the distance between the detector at the origin and the 
center of the ring, shell, and sphere at the points (0 0 

nr^Y^ " ^ '* HGURES 3 * 5 ' and 7 ' r «pec«iveiy. From 
HGURES 2-7. it can be appreciated that the detector must be. able to 
resolve three to four significant figures in order to reconstruct a map 
of 0.2 cm resolution with a field depth of 15 cm. 

Finite Detector Length 

The system function. of the Reconstruction Algorithm Section 
is the impulse response for a point detector. The following analysis 
will concern the impulse response for an NMR detector which has 
finite dimensions. Consider a detector array 401 as shown in FIGURE 
12 comprising multiple parallel planes 402 wherein each P !a nc has a 
plurality of antennae coils 403 wherein each has an area A The 
signal at any coil 403 is given by ,he integral of Eq. (29) over the area 
or the coi . This is the impulse response which is the system function 
which replaces h for a finite area detector. The Fourier Transform of 
this system function contains an argument of a product of the 
detector area A and the spatial frequency variables. Reconstruction 
could be performed as previously described in the Reconstruction 
Algorithm Section where this system function is substituted for the 
system function for a point detector. In the limit of zero area the 
measurement is that of a point sample. Thus, another approach is to 
use the lincanty of the superposition of magnetic flux to gain a hieher 
signal to noise ratio advantage by taking the difference of relatively 
large signals as opposed to performing the measurements with 
miniature antennas. In this case, the set of n coils each of area A of 
each plane comprise a grid of n blocks each of area a«A formed by 
the absence of common area overlap of one or more of ihe coils The 
signal of each element of the grid of area a can be solved by the set of 
linear difference equations of the signals of coils each of area A that 
corresponds to the noncommon area for the overlapping coils The 
solution of the signals due to the grid elements can be obtain using a 
computer by matrix inversion. 

The resulting values represent the average signal for each grid 
40 center locat.on^ The effect of this data processing operation on the 

spectrum can be modeled as a sample and hold, where the vohaecs at 
the centers of the grid elements arc sampled by multiplying by a 
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picket fence of delta functions spaced s units apart which are 
convolved w.th d, a square wave function in the x-direction and the 
l7^!° n ° f W ' dth S Uni,S Where ,he coordinates of Eq. (2) and 
FIGURE 9 arc used. In the frequency domain, this data processing 
5 operation causes the spectrum of the signal function s to be 

multiplied by D. the Fourier Transform, of the square wave function 
of width s units to form function S*. If this multiplication does not 
multiply S , the Founer Transform of the signal function, s, by zero for 
any frequency less than its bandwidth, then S can be recovered from 
10 S* by multiplying S* with the inverse of the Fourier transform of the 
sample and hold square wave function, a two dimensional sine 

ZT* , I y direcUons - This a " a 'ys* applies to all axes in 
which direction the detectors have finite length. Furthermore z- 

15 interval distances at wh.ch points discrete signals are recorded or by 
using multiple parallel plane detector arrays spaced a, the samplin* 
interval along the z-ax,s. However, if ,he signals are no, sampled a, 
discrete z-points. but each sample point is the integral resultant of the 
signal acquired continuously over a z-dispiaccment of q units which is 
20 much greater than the dimension of the detector ,n the z-direction 

then the corresponding sample and hold square wave has a width of □ 
units. " 



25 
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Altering The Dynamic Range Provided By The System Function 

The system function used in the Reconstruction Algorithm 
Section was the geometric system function, h, given as follows: 
h ~ 2z ~ * ~ v ' 

(x' + r + rf (65) 
This function is the impulse response at the detector array of the z- 
component of the magnetic field of a z-orientcd dipolc where the 
coordinates of Eq. (2) and FIGURE 9 arc used. The geometric system 
function convoluted with the advanced and delayed dipoles. which is 
represented by the function / (Eq. (34)) defined in the said section 
gives the function of the secondary flux when the tissue is magnetized 
35 with a constant magnetic field of unity. When the tissue is 

magnetized with a magnetic field whose function is not unity, the 
solution of Eq. (30) gives the magnetic moment of the voxel which is 
the product of die magnetic susceptibility and the external flux 
magnetizing the voxel. To obtain the magnetic susceptibility of the 
voxel, the solved magnetic moment is divided by the magnitude of 
the flux magnetizing the voxel. 
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Consider the NMR case where a quadratic magnetizing field 
gradient is applied along the y-axis in the z-dircction with the 
coordinates shown in FIGURES 10. 11, and 12 where the magnetizing 
field strength increases with distance from the detector array. A 
5 function of such a field is given as follows- 

B s -^*r:r (66) 
where a, and B 0 are constants and y o is the distance from the detector 
to the voxel. The system function in this case is given by the product 
of h and Eq. (66). And. the function of the secondary magnetic flux is 

10 given by the convolution of the function / (Eq. (34)) with the product 
of A and Eq. (66), where the advance of each dipole in the z-direction 
which appears in the function / (Eq. (34)) is substituted with the 
equivalent NMR coordinate, y. of Eq. (66). which reduces this part of 
the function to a constant. Thus, it follows that the solution of Eq (30) 

15 of the Reconstruction Algorithm Section is the magnetic moment of 

the voxel and that ihe magnetic susceptibility is given by dividing the 
said solution by the magnetizing flux. The purpose of applying a field 
gradient is to change the dynamic range of the signals of the 
secondary magnetic field. It was demonstrated in the Contrast and 

20 Limns of Resolution Section that the signal-to-noise ratio produced at 
the detector by the voxels determines the limits of resolution. The 
application of a gradient of magnetizing field which compensates for 
the drop-off of signal as a function of distance from the array 
improves the resolution relative to that which is possible in its 

25 absence. The gradient levels the signal as a function of distance. In 
effect, it prevents the contributions produced by voxels in close 
proximity to the array from dominating the detector signal such that 
contributions of distant voxels are lost in the relative error. Thus, the 
limits of resolution arc determined by the SNR of the detectors. 

30 In the absence of detector noise, the greater the variation of the 

signal as a function of position along the z-axis (y-axis in the NMR 
case), the greater the accuracy of the reconstruction. Noise 
considerations versus depth resolution arc weighed to determine the 
steepness of the magnetizing gradient if any which leads to an 

35 improvement in image quality. 

It will thus be seen that the invention efficiently attains the 
objects set forth above, among those made apparent from the 
preceding description. Since certain changes may made in the above 
constructions without departing from the scope of the invention, it is 

40 intended that all matter contained in the above description or shown 
in the accompanying drawings be interpreted as illustrative and not 
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in a limiting sense. In particular, optical methods of detection of 
magnetic resonance are within the scope of the present invention. 

It is also to be understood that the following claims are to cover 
5 all generic and specific features of the invention described herein, and 
all statements of the scope of the invention which, as a matter of 
language, might be said to fall therebetween. 
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APPENDIX I 

Integration of a RING of Dipolcs 

Derivation of the Field Produced by a Ring of Magnetic Dipoles 

The I -component of the magnetic field due to a dipoie or a loop of current of radius R 
with dipoie moment m = inR 1 at the origin is given as follows: 



B =j 



le-x-- v 



(M) 



The ? -component of the magnetic field at the position (x,y,z) due to a ring of dipoles of radius 
R with dipoie density -^L centered at the origin is 

((--vf Hv-y^U-^f 1 (L2) 

wherein (x',y\z') is a variable which corresponds to the position of each individual dipoie of the 
nng. 1 he relationsh ip bet ween Cartesian coordi nates and c ylindrical coordinates with t = 0 is 

= rtcos^ (i_3) 
y = /?sin# 

Substitution of Eq. (1.3) into Eq. (1.2) gives 

B ^UL i { 2zi -( x - Rc °s<P) 2 -{y-$tn<t>) 3 )j<i> 

Multiplying out terms gives 

fi w m V ( 2 ^ jy + ^ j cos->-2.r/gcos^) -(y- + ft'shV 0 -2y/?sm<p))</0 

2jr » ((r + R l cos" * - 2x/ccos0) + (y ! + ^ sin* 0 - 2yRs\^} + • = )** (L3) 
Substitution of Eq. (1.3) in the denominator and using the associative relationship and the 
trigonometric identity, cos' a + sur a - i m t he denominator gives 

B _ m 2 f { 2 z' ~{* } + & 'cos : <p ~2xRcos<t>)-( y 2 + r 2 m- <}>- 2yRsin<p))d^ 

z 2* I /V , i>l ^ 7 ~7Ts/7 ~ (1-6) 

o V +« -2/?(xcos# + ysm0)j 

The denominator can be approximated using the following relationship 
where 

^' ; +* ! (1.8) 
b = -2 tf(xcos # + y sin 0) {19) 

and 

/I = ~ 5/2 (1.10) 
This case gives the far field where « « 6 and/or r » ff. Using Eqs. (1.7-!. 10) gives 
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(r 2 + R 1 - 2/?(xcos £ + vsin )) 

I 5R 



Substitution of Eq. (Ml) into Eq, (1.6) gives 



7Tj7t(^cos^ + vsin< 



2n* 



I 5R , 



2*1 



In 

—I 



ft' cos 1 ft 2x/?cos0 



dtp 



5fo 3 cosft 5Rx l ys\n<t> SR^cos'f 
_ f 5£Vor »sin g lOx 2 R 2 cos' IQxyK 2 cos flsin g 



d</> 



5Ry ? .xcosft S/?/ sin ft SK*jcosftsin ' i 



2m 



SJfysin'ft jO^rvcos 0sin 0 lOtf V sin 2 <b 

" r ™ "T77T ~ " ■ 



The integration of Eq.(M2) gives 
4^ 2 2/cr 1 



2k 



xR l 



2*y 2 



(Ml) 



(M2) 



(M3) 



Collecting terms gives 
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Ring of Dipoles: 



B 



2l -X' - y 



sty-*/) 



(M4) 



a «„» J? (I . ,4 ;i Cmons,ra, « ^ ,hc ™g"«ic field of a ring of dipoles is not equal to that of 
» single d,pole at the ongin. The nng radius. R. appears in the denominator of each term The 

h! S5Tf ? . & ,S " r ° ™ C SCC ° nd and ,hird tcrms are additional P^rbaiions of 
the Held of a dipole at the ong.n whose magnitude is a function of the radius of the ring The 

SdlTlt ,h fi Crm : r VaniSh «"? When < he «*" ^ ^e nng is zero. A ring of dipolfs has a 
field that is cylmdncalty symmetrical. A shell, a cylinder, and a sphere of dipoles are the Z\ y 
other cases wh.ch have this symmetry. A cylinder is a linear combination of rings Thus the 
SJTSl ? ^ xt" dCm ° nS,fated b * sh0 --g *« it is different from that of a 
21 nf 5 1 3 SP ,; CrC - Thc prcScnt result ,ha ' t"e field of a dipolc is different for that of a 

All I r f ^ >nd APPENDIX IV. respectively, demonstrate that dipole field is un.que 
All other fields are a linear comb.nation of dipoles. Thus, the dipole is a basis element for the 
EEEEJ * Tr r SU t S " p,i r i,Uy ima8C - ThC -o'-iondepends on the e.temlt 
SJJ ° f ^ TI , ^' CS J d,ffCfS fr ° m lhat ° f 3 Sin 8' c di P° le al thc origin. The plot of the 
field ofanng of dipoles of rad.us R and magnetic moment m = J0< GW given by Eo (I 14) 

to th Tr" 0 " ° f R " " F1GURE 2 " ThC P~ ° f ' he thelg re a, vi 

to the detector is the po.nl (0.0.10). The plot of the field of a ring of dipoles of radius 

LzdZr d Z & T iCm0mCn u m = l °' Ccm> g^nbyEq-d l^asafunctionofthe 
distar.ee between the detector at the origin and the center of the ring at the points (00 

L , h Z"! Z ""I" $h0W " in FIGURE 3 Ff0m RGURES 2 and 3 it can be appreciated 

a 1 TJZT ™! * ab ' C 10 rc$0,VC * HrCC l ° f0Uf Si6nificant in ord " » ^onstmc, 

a map of 0.2 cm resolution with a field depth of 15 cm. 
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APPENDIX H 

Integration of a Sphericai Shell of Dipolcs 

Derivation of the Field Produced by a Shell of Magnetic Dipoles 

The z -component of the magnetic field due to a dipole or a loop of current of radius R 
with dipole moment m = inR 1 at the origin is given as follows: 

" 2z -x l -v- 

(Hl) 



B. =m 



The z -component of the magnetic field at the position (x.y.z) due to a shell of dipoles of radius 

R with dipole density - m ■ centered at the origin is 
4icR 2 

: = u ((x^)'*fr-yr*fc-^r <IL2) 

wherein (x\y'.z") is a variable which corresponds to the position of each individual dipole of the 
shell. The relationsh ip between C artesian coordinates and spherical coordinates is 

p -tJx 1 + v- + r r = v / } +y 2 +z 1 

x=psintf>cos0 x-tfsin#cos0 (113) 

y = psin<£sin0 y = RsmQsinO 

Substitution of Eq. (11.3) into Eq. (11.2) gives 

m yj(2(z-/?cos0) : -(x-/csm0cos£O 2 -^^ 

4rf oo ((x - Rsin #cos a) 2 + (y-tfsin *sin0) 1 cos 
Multiplying out terms gives 

B m 2 f'{ fa* -4zJtcQ$<p + 2R 2 cos 1 <l>-x 1 +2xRs\n<!>cQsG 
: * nR o 1 { x ~ ~ 2xflsin 0cos 0 + R 1 sin ! ^ cos 2 0 + y : - 2yrt sin 0 sin (9 

-ft* sin- ftcos* 8 - y 2 + 2y/?sin flsin 6 - R 1 sin ? flsin~ B)R 1 $\n<pd<t>d0 

+ R l sin : #sin- 9 r - 2zRcos<t> + tf* cos ; tf>) 5 ' 2 
Multiplying out the R 7 sin0 term gives 

B _ m l r | (2/?V sin 0 - 4zfl 3 sin 0cos 0 + 2 /?* cos" ^ sin <p - x : R* sin 0 + 2xR y sin" tfcos 0 
*" 4ffJ?' | J (x ; - 2x/?sin ^cosfi i ft' sin : £cos ; 0 + y- - 2 yrt sin 4 sin 9 

-R l sin' ftcos 2 9 - y'fl* sin 0 + 2vft 3 suv fl sinfl - ft 4 sin J ^sur 8\ " ^ 

-t-# 2 sin : ^sin : 0 + r - 2i/?cos0 + /? 2 cos : 
Substitution of Eq. (IL3) in the denominator gives 



(115) 
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4 ** 2 II (p* - 2xflsin ft cos a + R 1 sin 5 0 cos : 0 - 2y«sin ft sin 0 " 

-R l sin 3 ft cos 2 0 - v 2 /? 3 sin ft + 2yfl* sin* ft sin 0 - K* sin- ft shr 0\d<t>dd (U ' 7) 



+ 2 sin 1 ft sin 1 0 - 2zR cos 0 + R 1 cos 2 ft) 5 ' 2 

Using the associative relationship and using the trigonometric identity, cos : a +sin 2 a = I in 
the denominator gives 

B _ _m V r (2* V sin ft -4efl J s^ 

4 ^ II (p 2 + K 2 -2x/?sinftcos0 

sin* ft cos' 0 - y 2 R 2 sin ft + 2y/?' sin 2 ft sin 0 - fl 4 sin 3 ft s in 2 e]d<fx!e 8> 
~2 v/?sin ft sin 0 - 2z/ccosft) 5 ' 2 
The denominator can be approximated using the following relationship 

(a + b) m = a n + nba n ~ l (II.9) 

where 

<*=p 2 + ff ! (IMO) 
- -2xffsm £cos 0 - 2 vAsin ft sin 0 - 2z/?cosft (H. 1 1 ) 

and 

« = -5/2 (IL12) 
This case gives the far field where a « b and/or p» R. Using Eqs. (II.8-U. 12) gives 
(p : + # J -2xflsinftcos0-2yrtsinftsin0-2z/?cosft)" 5 ' 2 

B T~> — "^n*" *~ / ; 5 ^ M ?/i (^^^ft^^^ysinftsin(?-ncosft) < n ,3 > 

(p+R 1 ) (p + tf 2 ) r; 

Substitution of Eq. (11.13) into Eq. (IL8) gives . 

2ftVsinft JOtfVsinft, . 

(p 2 + R 2 )" 3 * [p 1 + R*) * Sin * C yS{n ^ $in ^ + zcos ^) 

4^sinftcos^ 2(W"! 4 sinftcosft / 

7^ — — (*sinftcos0 + vsinftsin 0 t £ cos ft) 
(P" + *J (p" + *) 

^cos^ftsin^ IOK J cos 2 ftsin0, . 

(p- + K J ) 5 '* + ' (p 2 4fl 2 )" 2 ( X5m ft cos fl + vsinftsin 0 + scosft) 

x 2 /? 2 sin ft 5x 2 /e 5 sinft, 

~71 — ^T^ 77 ~~ 71 7^7y ( x sin ft cos 0 + ysin ftsin 0 -f tcos ft) 

(p^/? 2 ) {p } +R l ) Y} 

2xR^sirrftcos0 lOxK 4 sin" ftcos 0 , 

Ia>^ n — (^sinftcos0 + vsinftsin0+ zcosft) 

(p" + /c) (p ) 

/? 4 sin j ftcos- 0 5K 5 sin- ftcos* 0 . 

{p 2 + R 2 f n [p 1 + R i ) in — Sm * C ° S * + VSin * S ' n ^ + Z C ° S ^ 

r^ 1 sin ft 5y 2 /? ) sinft / . 

7*1 — IToTT " 7"™; 7otT\ x sin ft cos 0 + ysin ft sin 0 + zcos ft) 

(p- + /r) (p- + tf 2 ) ' 



2*r 

B. = — — -y f f 
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2vK 3 sin 2 ftsin0 lOvtf 4 sin' ftsin 0 , . 

~~7^ pn Vi + — r~2 — ^T 7 ^ — (Jcsin^cos0+ysin^sinO+zcos^) 

IP + a J [P"+n I 

(-Tsinftcos0 + ysinftstn0 + zcosft) rf<W0 



fl 4 sin 1 ftsin' 0 3ft* sin* ftsin'' 0 



The integral of a sum is equal to the sum of the integrals. 



(II.I4) 



B. = — — -r{ || 



2iy S mfl JOtfVsin ft 



JJ 

0 0 

+jj 

o o 

-J J 

0 0 

0 0 
2* * 

-JJ 

0 0 
2ffX 

-JJ 

0 0 
0 0 

-JJ 



4 gfl 3 sin ft cos ft 20^ ' sin ft cos ft, 



2/? 4 cos : <frsinft lOtf'cos 2 ftsin ft. 

"T^ on*' 5 + "~T1 — I7n7r~(*sintfcos0 + y sin ft sin 0+z cos ft) 

IP + " ^ (P + A J 

x l R l smd 5 X 2 /*' sin ft, 1 

(p 3 + R*f r ~ + (p^T^j 777 * cos e + vsin ^ sin * + * cos ^ p*** 

2.rft J sin : <3cosfl 10*a* j sin ~ ftcosfl , 
(p 2 + fl 2 )*" + (p" + fl 2 )"* — ( tsm 0 c °s 0 + ysin ftsin 0 + z cos ft) Jc/ftrf0 

R 4 sin 3 ftcos-0 5K 5 sin 5 ftcos 2 0 



ttth — (xsin ft cos G + ysin ft sin 6 + ^ cos ft ) pft</0 



v 2 /? 2 sin ft 5 v 2 /t 3 sin ft, . ] 
/ , on sii + T~, — -77777 ft cos 0 + ysin ft sin B + z cos ft) Wft</0 
(P + * ) [P' + R) T 



ft)J*tf 



2y/? sin* ftsin g K) y/? 4 sin : ft sin 0 . 



, - - , , , T w «- . . v » v r I % o 1 1 i y/ 1 1 II %j . I 

(p' + R*)" 7 + (p- + R 2 y n — ( x5m * cos 9 + - vsin * sin 0 + * cos 0) p*/o 

/? 4 sin* ftsin- fl 5tf 5 sin 3 ftsin : 0, 1 ) 

(p* + * y ) vl *" "(p 3 + ^ YSm ^ cos0 ^ - VS2n ^ sin ^ + - c ^^)p^r (11.15) 



The integration of Eq. {II. 1 5) gives 
B = 



m 


Zn-R 1 


4nR 2 ' 





BkR* 



4ju : R 1 40/rrrt 4 



4nR > 4^R' Wm-R* 4nR' 

[P- 1 R'f ^P 1 * «T + >(P ! * *')'" ' 3(p J + Rf 1 



(H.16) 



Collecting terms gives 
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Shell of Dipoles: 



B. = 



m 
4;r/P 



40kR 4 



2z> 



-5* - .v 



(H.17) 



ftf ,, E . q . ( " ', 7) de ™ ns,ra, « ,hat ,hc '"agnetic field of a shell of dipoles is not equal to thai 

term Thl 31 C K rigin - ThC She " radiUS ' *' in « hc ^nLinatorof he fi7s t 

term. The first ,enn , S the magnetic field of a dipole at the origin only when the variab e 
correspond.ng ,o ,he radius of the shell. R. is KT0 . ^ secoi f d tc J is 

EST^^"*?? " * ° figin Wh0SC ™*™ udc is * function oft radius of 
me shell. I he second term vanishes only when the radius of the shell is yem Th„< ii~ a- , 

•s a basis element for the reconstniction of a magnetic su cen ibili^ ImZ tJ ^ 

S amfafaf r ieS ? f K rad ; US ^° 2Cm moment »,= 0< fa^ y 

Lq. II. n) as a funct.on of the distance between the detector at the origin and the center of 
shell at *c po,nts ((M>. z) is shown in FIGURE 5. From FIGURES S 5 iic^T 
appreaatcd tha, the detector must be able to resolve ,htee to four significant figures in order to 
reconstruct a map of 0.2 cm resolution w,,h a field depth of 15 cm 
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APPENDIX III 

Proof that the Field Produced by a Shell of Magnetic Dipoies is Different from that of a Single 

Dipole 

Consider Eq. (H.5). 

B _ m Vj (2z ; -4jj?cos0-h2/? :? co$ 2 ^^jr ; ^2x/?sin^cosg 
* KR * yJP ~~ 2.rKsin 0co$0 + R 2 sin* tfcos 1 0 + y l ~ 2yRs\n<j>s\n6 

-R 2 sin : (pcos'O-y 2 + 2yffsin<Mn 0 - g sur ftsin 2 g)*' sin fr/fr/fl lHLI> 
+ /? 2 sin 2 £sin 2 0 + r - 2z#cos^ + R* cos 2 $f l ~~ ~~ 
The conditions for this integral to equal that of the field of a dipole at the origin are 
^z*cos 0 + 2tf 2 cos-> + 2*/?m 

and < IU2 > 

-2xRsm 0 cos 6 + /?- sin : 0cos ; 0 - 2yKsin 0sin 0 + tfW>sm J 0 - 2zffcos0 + # 2 cos : 0 = 0 

Thus, Eq. (f IL2) must equal Eq. (III.3). 

-4*/?cos0 + 2* 1 cos- 0 + 2*/?sin 0cos0 - tf- sin : 0cos 2 0 + 2 vtfsin 0sin 0 - /?' sin 2 fsin 1 0 

« -2*i?sin 0cos0 + /?-sin>cos ? 0 - 2yrtsin 4>sin 0 4- /<' sin 2 ^sin 2 0 - l~JRco$$ + /c*cos 2 0 

T , r (IH.4) 
Therefore, 

~2zRco$<p+ R l cos 2 0 - 2 tf 2 sin : 0cos J 0 + 4 x sin 0cos 0 + 4yRsm 0$in 0 - 2K 2 sin ? 0sin 2 0 = 0 

., • t (ni-5) 

Using the associative relationship and the trigonometric identity, cos : a + sin" a - 1 as well as 
factoring out an R in Eq. (HI.5) gives 

-2zflcos * + ! cos 2 0-2/?' sin' <p + 4 * sin tfcos 0 + 4y/?sin 0sin 0 = 0 (IH.6) 
-2zcos0 +3tfcos>-2/? + 4.vsin0cos0 + 4ysin0sin0 = 0 (HI.7) 

For x * 0. y * 0, z * 0. Eq. (IH.7) is true only if R = 0 which proves that the field of a shell of 
dipolcs centered on the origin is different from that of a single dipole at the origin. 

Eq. (HI.5) can not be integrated in closed form; however, the integral can be 
approximated so that the deiector tolerances for a given image resolution can be determined. 
The denominator of Eq. (HI.5) can be approximated using the following relationship 

(a+b) m ~a*+nba~ l (m 8) 

where 

a=p 2 + R> (ni9) 
b ss -2*rtsin tf>cos 0 - lyRsin 0sin 0 - 2zRcos<p (INTO) 

and 

^ x n = -5/2 (mn) 
The result given by Eq. (U.17) is 
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4/rK 2 



UnR 



2v 



(p- + R) J 3 



3 [(p' + R>) m 



(III. 12) 



Higher order terms of the approximation given bv Eo mi «\ . • - , 
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APPENDIX IV 



Integration of a Sphere of Dipoles 

Derivation of the Field Produced by a Sphere of Magnetic Dipoles 

The z <omponent of the magnetic field due to a sphere of dipoles is derived using the 
equation for a shell of dipoles with the substitution of the dipole density of the sphere 3m } 

where the radius of the sphere is R. The field of a shell of dipoles is given approximate!^ by 
Eq. (II. 1 7) of the Derivation of the Field Produced by a Shell of Magnetic Dipoles Section 



B = 



4/rr 



4*r : 



40ro-' 



2z> 



* l -y' 



) 



(IV. I) 



where the radius of the shell is r (The parameter r replaces the constant 'r of Eq. (II. 1 7) as 
the radius of the shell.). For small r. the first term of Eq. (IV. I ) dominates and is used to 
approximately calculate the z -component of the magnetic field due to a sphere of z -oriented 
dipoles as follows: 

A sphere of dipoles is equivalent to the integral over concentric shells of dipoles each 
of radius r where 0<; r < R. Thus, the i -component of the magnetic field due to a sphere of 
dipoles is the integral of the field of the shells given by Eq. (IV. I). 



B. = 



T 4*J 



Let 



2z l 



r*dr 



= tan 0; dr = psec 2 Odd 



and using the trigonometric identity 

1 + tan 2 9 = scc ! 0 
Eq. (IV.2) becomes 

B .-=(2z'-* ! ->'-) T ~| 



Substitution in of 

sinO 



4k tan 1 6 

3 



sec' OdO 



and 



cos0 



1 



sec 9 = 

CQSO 

into Eq. (IV.5) gives 



m 



Let 



3 



f *n » ^ 

—smOcosOdO 

i P~ 



(IV.2) 

(IV.3) 
(IV.4) 
(IV.5) 

(IV.6) 
(IV.7) 

(IV.8) 
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u = sin 0; du - cos 0d$ 
Substitution of Eq. <IV.9) into Eq. (IV.8) followed by integration gives 

3 

Substitution of Eq. (IV.9) into Eq. (IV. 10) gives 
B.=(2? , -jr-.v s ) m 4 * sin '* 



From Eq. (IV.3) 



3 



- = tanfl; sin9 = - — £- 
P 



Substitution of Eq. (IV.12) into Eq. (IVJ |) gives 



fly/ , oM'Ho 
1+1 



Evaluation at the limits of the integral gives 



B : =(2r-.T i -y > ) 



m 



Multiplication of Eq. (IV. 14) by 1 



Hffl 
US) 

K-en 



gives 



Multiplying out the p 5 term in the denominator gives 



UV.9) 
(IV. 10) 

(IV 11) 

(IV.I2) 
(IV.I3) 
(IV.I4) 



(IV. 15) 



SUBSTITUTE SHEET (RULE 26) 



90 



Sphere of Dipolcs: 



(IV. 16) 



.h,, Eq (I , V 'i 6) derno " s,nucs ,hat ,he m^netic field of a sphere of dipoles is not equal to 

EL The£ Slsl origin - T hc M sheU radius - *• appcars in ,he o L 

term. The first term is the magnetic field of a dipole at the origin only when the variable 

penurbat,on of the field of a d.polc at the ongin whose magnitude is a function of the radius of 
the sphere rite second term vanishes only when the radius of the sphere is zero , TW ,Z 
dtp* , a bas.s element for the reconstruction of a magnetic «Z£ '£1 T h ' 

FICJIPP/; Tk* „ rT gI c " D y^UV.l6)asarunctionof radius /? isgivenin 

(0 0 of The m r CCm r ° f *" Sf>hCfC rda " Ve t0 ,hc dc<cc «» is the poin, 

(U.U.10). The plot of the field of a sphere of dipoles of radius R - 0 2 r m a „,t m „ . 
moment „, = I0 4 Gcm J R iv en bv Eo iIV « , r ■ I H " magnetic 

fi E „L bC , aPPreC ' a,Cd ** < he dc '«'°' ™* be ab.e to resolve three toTur 

S.gmf.cant figures ,n order to reconstruct a map of 0.2 em resolution with a field depth of 1 5 
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APPENDIX V 
Fourier Transform of the System Function 
The system function. h( P ,<p, z ), j„ cylindrical coordinates is 

(x+v+r)"-' |p- + ^f (V.I) 

The spacetime Fourier transform in three dimensions in cylindrical coordinates //(* <i> * ) is 
given (9} as follows: * rn V <p '*->- 1S 

- it- 

H(k, .<».*.) - J J J „(p, z )exp(-/24^cos(d, - «>) + * ;? j)p</p^ ( v 2) 

Wiih circular symmetry [9} 

H{k,X) = 2nj jh(p, z )J AkpP)e -^ pdpdz {v 
The Fourier transform°o7the systen, function is given by the substitution of Eq. (V. , ) int0 Eq. 

h <.h$3f j a**w*-*' 1 * (V4) 

Consider the integral of Eq. (V.4) wi.h respect to dp only. Factorization of *(p.*. c ) gives 

Consider the definite integral 

r- r^UflfJrfr _ gV^jj^Jazj 

>° [r+z l \"> " 2T|m+1] (V.6) 
and the modified Bcsse! function of the third kind relationship 

The first factor of Eq. ( V.5) is the same form as Eq. (V.6) with v = 0; u = f ,hus. 
same fonn as Eq. (V.6) usmg the Besse. function of the first kind recurrence relationsh.p 

'-"♦'~M-f''M (V . 9) 

Consider the second factor of the integral of Eq. (V.5). 

-"t&y^w «v,o, 

Eq.(V9)with U=| is 

x ,l 1 (V II) 
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Let 

Substitution of Eq. (V.I3) into Eq. (V.I2) is 

Substitution of Eq. (V.IO) into Eq. (V.I 4) is 
f P 



(V.I2) 
(V.13) 
(V.!4) 





2 







(V.I 5) 



The first factor of the right hand side of Eq. (V. 15) is the same form as Eq. (V.6) with 

v=*l; « = thus, 
2 

where K_ v3 [k p z) » K m [k p z] (Eq.(V.7)). The second factor of the right hand side of Eq. (VA5) 

is the same form as Eq. (V.6) with v = 2; u = - thus 

2* 

2 ^^T7r J ^ pvp =^^ tv.17) 

Combining the pans of the integration with respect to rfp of Eq. (V.4) by adding Eq (V 8) Ea 
(V.16),andEq.(V.17)givcs * H ' 

"{2" J |nz , "jr'' 3 [2" l W" J Jt" 1 1 



J. 



The modified Besscl function of the third kind formulae is 

L2xJ m!r(n + l-m| 

Substitution of Eq. (V. 13) into Eq. (V.I9) with u = 1 is 



2^ 



Substitution of Eq. (V.13) into Eq.(V. 19) with v = 0 is 
Substitution of Eq. (V.20) and Eq. (V.2I) into Eq. (V.I8) is 



(V.18) 
(V.I9) 

(V.20) 

(V.21) 
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H5/2J 



2k f z 



H5/2J + (2" l )r(5/2]|2^J ' 



H5/2)2 



(V.22) 



H5/2J 



r{5/2)2 



Collecting terms gives 

With r(3] = 2a<id r(5/2] = 3/4s" t ( Eq.{V.24)is 



A- 7T )/Z 



Integration of Eq. ( V.29) with respect to d Z gives 



4tl* 



*ives 



m of Eq. (\ 



Multipiicalion of Eq. (V.3 1 ) by 
1 = 



Arck. 



-jk z + k p 



k: +t 



(V.23) 
(V.24) 

(V.25) 

(V.26) 

(V.27) 
(V.28) 
(V.29) 

(V.30) 
(V.3I) 

(V.32) 

(V.33) 



The system function (Eq. (V.I)) is an even function; thus, the spacclime Fourier transform in 

). is given by taking the real pan of Eq. 

(V.33) (8). 

' • v+v < v - 34 ) 

The spacetime Fourier transform in three dimensions in Cartesian coordinates //(A A) is 



I*/ 



(V.35) 
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where «he relationship between «he wavcnumbcrs and the spatial Cartesian coordinates is ; 



follows: t—w «-« <-iwrainaics is as 

A, x (V.36) 

2£ 2£ 

A. " v (V.37) 
2*_2* 

A. • r (V.38) 
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APPENDIX VI 
Derivation of S^HF^Ufc) from Eq. (44). 

. . Ink 7 ~ r> w >' 2 

where the Fourier transform of = i/<?r z > 0 a/«/u(*) = I /or z < 0 (8) is 

= + (VI. 2) 

and where // is given by Eq. (41) and F is given by Eq. (38). The convolution integral of the 
second term of Eq. (VI. I) is 

V" ***** -/I*. '.•*.» J f 3^*. -.^ 

(VI.3) 

Collecting exponential terms gives 



" Meeting exponential terms gives 



S=- 

(VL4) 

Expansion of the denominator of the second rerm gives 

(VI.5) 

The factorization of the denominator of the second term, k] - 2k, k. + + A?, using the 
quadratic formula is 



2^±[4^-4[^+A:/]p 
_ 

2A.±[4^-4^-4^j ? 



(VI.6) 



2 

K z =k.±jk, 
Substitution of Eq. (VI.6) into Eq. (VI.5) is 

(VJ.7) 

The expansion of Eq. ( V1.7) by the method of partial fractions is 
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2 I X. e'^""^-\ 



B= 1 



(VI.8) 



(VJ. 12) 



J L [k z -* 5 IT. - y* p 

The factors A. B, C of Eq. (VI.8) arc determined as follows: 

Let k : = 0. Then, Eq. (VI.9) gives 

A [ k P + 1 (VI. 10) 

Thus, 

Lei K : =k. + jk F . Then Eq. (VI.9) is 

+ C[k z +jk p lk z +jk p -k z +jk fi }=\ 
The term in B vanishes. Substitution of Eq. (VI. 1 1) into Eq. (VI. 12) gives 
[5^J][*- ? + W, -2*/ -2jfc*, + *; +k? ) + C(2;V, -2*/]=i ( vi.| 3 ) 
where the term in /\ also vanishes. Solving for C gives 

C= R?^S] (VII4) 

Let K\ = k. - jk p Then, Eq. (VI.9) is 

The icrm in C vanishes. Substitution of Eq. (VI. 1 1) into Eq. (VI. 1 5) gives 
fff^jl*' - 2 A*, " ^ - 2* : ' + 2yU p + + fi[-2,* * p - 2^)= 1 (VI. 16) 

where the term in A also vanishes. Solving for B gives 



(VI. 1 5) 



Substitution of A, # f and C into the convolution integral (Eq. (VI.8) gives: 
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Ink 



.i - 

+— M j S . 

The first convolution integral is of the form given by Mc Scibert J8J. 
— <=> sgnz n \ sgrtz, ~ -I where z h < 0 

Thus, 

[77 777 j ~~ ~ ~?~7 o w hcre <0 

|% +*:]/. J*, \ k P + k z\ 

A change of variable in the second convolution integral given by letting 



tn 



JS 



J ~ l I"*-} where z n <0 



-e~ M:J where r <0 



A change of variable in the third convolution integral given by letting 



tn 



is 



Ink: 
p 



K' + 2;^)A 4-2*J whcrc *- <0 



(VI. 1 8) 



(VI. 19) 

(VI.20) 
(VI.2I) 
(VI.22) 

(Vl.23) 

(VI.24) 

(VI.25) 
(VI.26) 
(VI.27) 

(VI.28) 
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2zk;e 



J ~r—dk. 



|-2*;+2jU,]A 2 A 



(VI.29) 



f e' M:J where - <0 



2|*;-A*p| — (VI.30) 

Combining Eqs. (VI.20). (VL25). and (VI.30) gives the convolution of Eq. (VM) 
I • "'• 4 '-.n ♦«,/? r , » ^ ' 

5 4 wf+ i i i 



"fx 



W,/2 



(VI.3!) 
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APPENDIX VII 

Derivation of chc Inverse Transform of Eq. (57) to Give Inverse Transform 1 , Eq. (58). 



|«:.*0 [ k p \ I 

2* J j C '^p|A// ; * 



-'. >t-n 



+2n j } e'^'yjj^j^e^ 

— 0 



where : <0 



Consider the definite integral: 

o |<i + f>j ;r 2 

In the case of the first integral of Eq. ( VH.2). the parameters of Eq. (V1I.3) 

v = 0 

' = *, 
b = p 



arc 



Thus. 



where 



Consider the following relationship of the modified Bcssel function of the third kind 

«tM..* 4'" ; h'-'" K P' : ; 

In the case that v - 0, Eq. (VII.7) is 
Thus, 



(VIM) 



(VII.2) 
(VH.3) 

(VII.4) 

(VII.5) 

(VH.6) 

(VH.7) 
(VII.8) 
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With a change of variable as follows: ' V "" 
First Integral: Second Integral: 

and an expansion of the complex exponential factor, Eq. (Vn 10) becomes 

= ^j|(cos^.J +/ sin[t pJ J + cos[vJ-ysin(^^])^(A p pj^ (V „ |2) 

The transform is given by Bateman {12) 

_, (VH.13) 

Thus. 

The Fourier transform relationship between a function f\x\ and its Fourier transform e\ v] 
given by McC. Sieben (8) is 

*A*\~--.j ji m (VIIl5) 

Consider the following Inverse Fourier transform integral from Eq. (VII 2) 

J V'-'^KpK (VH.I6) 

0 

Eq. (VIU6) is equivalent to the following Inverse Fourier transform 

^'W*4 (VII.I7) 
From Eq. (VII. 15) with 

^"R+i^f (VII,8) 

Eq. (VII. 16) is 
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rcr/us98tt<i22o 



-i 



(VII. 19) 



Tak.ng Ihe imagtaiy pin bcca.se Eq. (VII. 19) has odd symmeicy gives 

jv'-'^K'p^r-^O (V1U0) 



Thus, 



(VII.2 1) 



The Inverse Fourier .ransform with respec. ,o * s is given by MeC. Sieben |8|. In ,hc ease of 
the first two terms of Eq. (VII.2), 

In the case of the third and fourth terms of Eq. (VII. 2), 
J jky^ : dk : = S[z) 

second terms of Eq. (VII.2) and comb.n.ng transforms from Eqs. (VII. 14) (VH.22) and 
(VI .24) ,n the case of the third and fourth terms of Eq (VII.2) gives 

Hi 2 I £«. «•* 



(VII.22) 



(VII23) 



(VH.24) 



= r 



8 ^ *«ii.«j^>'" 



2. I I *. e-'l'— ' -I 



e''l-' +e -^.l + y £,| / M:.|_ f - M; .,| 



where < 0 

Substitution of Eq. (37) into Eq. (VI1.25) gives 



-1 



(VH.25) 
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♦/./> ./,/? 



2 



+ 2/r£[zj 



-I 



where * <0 



Taking (he derivative given by ihe doublet function gives 

i t 1 x. i^M e -^ + e .*i + jLi e '*s _ e -, 

? ..-f : /2 ,.,-/,/> ' k t 



(VII.26) 



" 7 Z S S Jf-.ji.^.*( jr -"i*.y-«.*.z-n3Jtj® 
Inverse Transform I (VJI.27) 



S[i) 
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What is claimed is: 



I • A system for providing a multi-dimensional magnetic susceptibility image of a 
body, comprising: 

a first radiation source for applying a first radiation field having a magnetic 
field component to magnetize the body, 

a first detector for measuring said magnetic field in absence of the body in a 
volume to be occupied by the body, 

means for applying a second radiation field to the body when disposed in 
the magnetic field to produce a third radiation field emanating from the body, 

a second detector producing output signals in response to said third 
radiation field, and 

a processor coupled to said second detector for receiving said output signals 
to produce a multi-dimensional magnetic susceptibility image of the body. 

2. The system of claim I, wherein said radiation source comprises a magnet for 
producing a DC magnetic field. 

3. The system of claim 2, wherein said magnet comprises a superconducting magnet. 

4. The system of claim I , wherein said first detector comprises a magnetometer. 

5. The system of claim I , wherein said means for applying a second radiation field 
comprises an RF generator for applying an Rf pulse to the body. 

6. The system of claim I, further comprising an analog to digital convener coupled to 
said second detector for digitizing said output signals. 

7- The system of claim I, wherein die body comprises at least one type of nuclei 
having non-zero spin, and wherein said second radiation field polarizes said nuclei and 
said polarized nuclei precess to produce said third radiation field. 

8. The system of claim 7, wherein said at least one type of nuclei comprises protons. 

9. The system of claim I, wherein said processor includes a Fourier transform stage 
for obtaining Fourier spectra of said output signals of said second detector. 
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10. The system of claim 1, further comprising a plurality of detectors for detecting said 
third radiation field, said plurality of detectors being disposed in a plane with respect to 
each other and with respect to said second detector, said plurality of detectors being 
configured to move in a direction perpendicular to the plane to sample said third radiation 
field over a three-dimensional volume external to the body. 

11. The system of claim I, further comprising a plurality of detectors forming a three- 
dimensional array to sample said third radiation field over a three-dimensional volume 
external to the body. 

1 2. The system of claim 10, wherein said processor comprises a Fourier transform 
stage for obtaining frequencies and intensities of Fourier components of said sampled third 
radiation, said frequency components indicating the magnetic susceptibility corresponding 
to voxels of the body and said intensity components indicating positions of said voxels. 

1 3. The system of claim t. wherein said third radiation field emanating from a voxel of 
the body in said magnetic field is a near field of a dipoic, wherein the near field has a 
maximum amplitude given by the following formula: 

(jr + jr +i y n 

wherein x.y, and z are Cartesian coordinates, and rn 2 is a bulk magnetization of the 
voxel along the z-coordinate, and B* is the maximum amplitude of the near field. 

14. The system of claim I, wherein said processor includes a Fourier transform 
processor for obtaining. Larmor frequencies corresponding to voxels of the body from 
which said third radiation field emanates. 

1 5. The system of claim 1 , further comprising a display element for generating an 
anatomical image of the body from the multi-dimensional susceptibility image. 

16. The system of claim 14, wherein said display clement generates said anatomical 
images based on selected physiological parameters, said anatomical images being 
substantially free of motion artifacts. 
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17. System for providing a multi-dimensional image representation of spatial 
variations of magnetic susceptibility of a body including a paramagnetic or a diamagnetic 
substance, comprising: 

a first generator for generating an excitation field in the body to polarize 
5 paramagnetic or diamagnetic substances within the body such that said polarization 

creates a local magnetic field in voxels of the body according to the magnetic 
susceptibility of each voxel, 

a second generator for generating a primary radiation field in the body to 
induce a secondary radiation field to emanate from the body, each voxel of the 
10 body contributing to said secondary radiation field at a frequency determined by 

the magnetic susceptibility of said voxel, 

a detector positioned to receive one of said fields for producing an output 
signal in response to said secondary radiation field, and 

a processor for receiving and processing the output signal of said detector 
15 to create the multi-dimensional image. 

18. System for providing a multi-dimensional image representation of spatial 
variations of magnetic susceptibility of a body having at least one type of nuclei with non- 
zero spin, and a paramagnetic or a diamagnetic substance, comprising: 

20 a first g cncra '<>f for generating an excitation field in the body to polarize the 

paramagnetic or the diamagnetic substances such that said polarization creates a 
local excitation field in the body, said local excitation field generating a shift in the 
Larmor frequency of nuclei in the body, 

a second generator for generating a primary radio frequency field in the 
25 body to create a transverse magnetization of the nuclei, 

a detector for producing output signals in response to detecting a 
component of a secondary radio frequency field generated by the polarized nuclei, 
and 

a construction processor for receiving the output signals of said detector 
30 'o create the multi-dimensional image. 

1 9. System for generating a three-dimensional magnetic susceptibility map of a body, 
comprising: 

means for producing a magnetic field for magnetizing the body; 
35 means for measuring the magnetic field in the absence of the body in a 

volume to be occupied by the body; 
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a radiation generator and a transmitter for producing a primary radiation 
field for exciting the body: said excited body producing a secondary radiation Held 
external to the body; 

a planar array of detectors movable in a direction perpendicular to the plane 
of the detectors for producing output signals by spatially sampling a selected 
component of said secondary radiation field over a three-dimensional volume at 
least at the Nyquist rate; 

an analog to digital converter for digitizing said output signals; 
a Fourier transform stage for converting the output signals of the detectors 
acquired over a period of time into frequency spectra; and 

a processor for converting the frequency spectra into a map of magnetic 
moments in the volume and for converting the magnetic moment map into a three- 
dimensional magnetic susceptibility map. 

1 5 20, The system of claim 1 9. further comprising an image processor for displaying said 
magnetic susceptibility map, said processor being configured to display said map from any 
perspective as a two-dimensional or a three-dimensional image. 

2 1 . The system of claim 1 9, further comprising coils for providing a magnetic field 
20 gradient to alter intensity of the said secondary radiation field when said radiation 

generator provides Tj and T? radio frequency pulse sequences. 

22. The system of claim !9 X wherein said radiation generator comprises a radio 
frequency generator providing a radio frequency field to cause precession of protons in 

25 said body to produce said secondary radiation field. 

23. The system of claim 19. wherein said means for measuring the magnetic field 
comprises a magnetometer that employs NMR frequency of protons in water to detect 
magnitude of said magnetic field in the absence of the body in the volume to be occupied 

30 by the body. 

24. The system of claim 1 9. wherein the body comprises substantially non- 
ferromagnetic material. 

35 25. The system of claim 19. wherein the radiation generator includes an amplifier for 
increasing magnitude of said primary field in portions of the body distal to said array of 
detectors. 



SUBSTITUTE SHEET (RULE 26) 



WO 99/26078 



107 



PCT/US98/24220 



26. The system of claim 1 9. wherein said means for producing a magnetic Held 
produces a magnetic field varying as a quadratic function of distance from said detectors in 
a direction perpendicular to the plane of the detectors. 

5 

27. The system of claim 19. wherein said magnetic field has a component B z in a 
voxel along a Cartesian coordinate perpendicular to said array of detectors, said 
component B z being given by: 

B 2 = B Q [ a 2+ y 2 m |3/2 whcre a(> and Bo are conslants and h {hc dislance fr0m 
10 the plane o f the detectors to the voxe I 



28. The system of claim 19. wherein said magnetic field :z substantially 



constant. 



29. The system of claim 19. wherein said magnetic field is substantially confined to a 
1 5 volume of the body to be imaged. 

30. The system of claim 19. wherein said means for producing a magnetic field 
produces a magnetic field having a magnitude that varies within the body, said varying 
field providing a range of Larmor frequencies. 

20 

31. The system of claim 19. wherein said array of detectors comprises an array of 
antennae responsive io a selected component of said secondary radiation field. 

32. A method for providing a multi-dimensional magnetic susceptibility image of a 
25 body, comprising the steps of 

measuring a magnetic field component in a volume to be occupied by the body, 
applying a first radiation field having the magnetic field component to the body to 
magnetize the body, 

applying a second radiation field to die body when disposed in the magnetic field 
30 to produce a third radiation field emanating from the body, 

detecting the third radiation field emanating from the body. and. 
producing the magnetic susceptibility image of the body from said detected third 
radiation field. 



35 



33. The method of claim 32. wherein the step of detecting comprises the step of 
positioning an array of detectors movable in a direction perpendicular to the body to 
sample said third radiation field over a three-dimensional volume. 
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34. The method of claim 32. wherein said step of applying a first radiation field 
comprises the step of applying a DC magnetic field to the body. 

35. The method of claim 32. wherein the step of applying a second radiation field 
comprises the step of applying an RF pulse to the body. 



36. The method of claim 33. wherein said step of producing the magnetic susceptibility 
image further comprises the steps of 

10 determining the frequency components of voxels of said third radiation field. 

determining the magnetic susceptibility of said voxels at each said frequency 
component, and 

obtaining the positions of said voxels by utilizing intensity variation of said 
frequency components over the sampled three-dimensional volume. 

15 

37. A method for providing a multi-dimensional image representation of the spatial 
variations of magnetic susceptibility of a body having a paramagnetic or a diamagnetic 
substance, comprising the steps of 

generating an excitation field in the body to polarize the paramagnetic or the 
20 diamagneiic substance to create local magnetic fields in voxels in the body according to 
the magnetic susceptibility of each voxel, 

generating a primary radiation field in the body to induce a secondary radiation 
field from the body, each voxel of the body contributing to said secondary radiation field 
at a frequency determined by the magnetic susceptibility of said voxel. 
25 generating an output signal in response to said secondary radiation field and 

indicative of said radiation field, and 

creating the multi-dimensional susceptibility image from said output signal. 

38. A method for generating a three-dimensional magnetic susceptibility map of a 
30 body, comprising the steps of 

producing a magnetic field for magnetizing the body. 

measuring the magnetic field in the absence of the body in a volume to be occupied 
by the body. 

producing a primary radiation field for exciting the body, said excited body 
35 producing a secondary radiation field external to the body. 

disposing a planar array of detectors external to said body, said planar array of 
detectors being movable in a direction perpendicular to the plane of the detectors for 
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producing output signals by spatially sampling a selected component of said secondary 
radiation Held over a three-dimensional volume at least at a Nyquist frequency, and 
producing the magnetic susceptibility image of the body from said sampled 
selected component of said secondary radiation field. 

5 

39. A method for generating a three-dimensional magnetic susceptibility map of a 
body, comprising the steps of 

producing a magnetic field for magnetizing the body. 

measuring the magnetic field in the absence of the body in a volume to be occupied 
10 by the body. 

producing a primary radiation field for exciting the body, said excited body 
producing a secondary radiation field external to the body. 

disposing an array of detectors external to said body, said array of detectors being 
movable in a direction perpendicular to the detectors for producing output signals by 
15 spatially sampling a selected component of said secondary radiation field overa three- 
dimensional volume at least at a Nyquist rate. 

digitizing said output signals. 

convening said output signals acquired over a period of time into frequency 
spectra. 

converting said frequency spectra into a map of magnetic moments in the volume 
and convening the magnetic moment map into a three-dimensional magnetic susceptibility 
map of the body. 



20 



40. A method for producing a multi-dimensional image representation of the spatial 
25 variations of the magnetic susceptibility of a body having a paramagnetic or a diamagnctic 
substance, comprising the steps of 

exciting said paramagnetic or diamagnetic substance by applying an excitation 
field to said body. 

detecting a radiation field originating in said body from said excited paramagnetic 
30 or diamagnetic substances, wherein frequency components of said radiation field provide 
information regarding the variations of the magnetic susceptibility within said body. 

constructing the multi-dimensional susceptibility image according to the 
information provided by the frequency components of said radiation field. 

3a 41. A method for providing a multi-dimensional image representation of spatial 
variations of the magnetic susceptibility of a body, comprising the steps of: 
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providing a primary magnetic Held to magnetize said body, 
measuring said primary magnetic Held in the absence of said body in a volume to 
be occupied by said body. 

applying a primary radio frequency field to said body to rotate at least one type of 
magnetic moments in said body by 90 degrees such that said rotated moments point in a 
direction transverse to said primary magnetic field. 

detecting a secondary radio frequency field produced by a free induction decay of 
said rotated magneiic moments by a plurality 0 r detectors arranged in a plane and movable 
in a direction perpendicular lo said plane. 

obtaining Fourier transform of said secondary radio frequency to derive Larmor 
frequencies of said rotated magnetic moments. 

determining the magnetic moment corresponding to each Larmor frequency. 

determining maximum intensity of each Larmor frequency at positions of said 
detectors, and 

constructing the multi-dimensional susceptibility image according to said derived 
Larmor frequencies and according to spatial variations of said secondary radio frequency 
field. 



42. The method of claim 4 1 . wherein said constructing step comprises the step of 
20 employing a reconstruction algorithm. 

43 . A method for determining a correspondence between detected voltages in an array 
of detectors lo a radiation distributed within a body, comprising the steps of: 

applying a magnetic field to the body. 
25 applying a first radiation field to the body to cause the body to cause a radiation 

distribution within the body, 

detecting the radiation distribution within the body with an array of delectors. and 
determining the correspondence between voltages of the array and the radiation 
distribution such a set of detected voltages corresponds to a unique radiation distribution 
30 within the body. 

44. The method of claim 43. further comprising the step of constructing a magnetic 
susceptibility map of the body in response to one or more voltages detected by said array 
of detectors. 

35 

45. The method of claim 43. further comprising the step of correlating detected 
voltages with a particular frequency of a magnetic field emanating from the body and 
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induced by said applied field, and determining the spatial locations within the body 
contributing to each detected frequency. 

46. The method of claim 45. further comprising the step of constructing a magnetic 
5 susceptibility map of the body as a function of said detected voltages. 



SUBSTITUTE SHEET (RULE 26) 



WO 99/26078 



PCT/US98/74220 



1/15 



110 



124 



122 



126a 



RECOMSIRUCIIOM 
PROCESSOR 



ff 

PROCESSOR 



SPAIIAL 
I OCA I OR 



-126 



-126b 



A/D 



A. 



DRIVE 
MECHANISM 



OlSPLAt 
PROCESSOR 



-I2S 



112- 




crarcraacjao 
o a a cj a a a o 
a cf a a o cr a a 



118- 



716a V/7 



5 



o 




-> z 



M ACME 10 
METER 



-116 



RE 
CEMERAIOR 



FIG. la 



SUBSTITUTE SHEET (RULE 26) 



PCTVUS98/I4220 



2/15 



228 




-> z 



FIG. lb 



SUBSTITUTE SHEET (RULE 26) 



WO 99/26078 



PCT/US98/24220 




FIGURE lc 



WO 99/2*07$ 



PCT/US98/24220 



4/15 




FIGURE 2 



WO 99/76078 



PCI7US98/24220 



5/15 




Distance (cm) 



FIGURE 3 



WO 99/26078 



PCT/US98/24220 



6/15 




FIGURE 4 



WO 99/26078 



PCT/OS98/24220 



7/15 




FIGURE 5 



WO 99/26078 



PCT/US98/24220 



6/15 



Radius (cm) 




FIGURE 6 



WO 99/26078 



PCT/US9S/24220 



9/15 




FIGURE 7 



WO 99/26078 



PCTAJS9 8/24220 



10/15 



25- 



fit to 

H£ASU8fM(AfT 

AA/O 
£8R0X SfCi/Al 



StCNAl 



-22 



X 

So 

COMTAOL 




24 
V 



PXOCLSSiNG 
AAfO 

CO/fTXOL CIRCVtTS 



23 



01 SPLAY 



FIGURE 8 



WO 99/26078 



PCT7US98/24220 



11/15 



Z 




FIGURE 9 



WO 99/26078 



PCT/US98/24Z20 



12/15 



y 




z 



X 



FIGURE 10 



PCT/US98/24220 



13/15 



\ 



203 




14* 



Voxel Position Relative to the 
Detector ot the Origin of the 
Sample Space is (0,0,-a) 



FIGURE 11 



WO 99/26078 



PCT/US98/24220 



14/15 




403 



404 



Voxel Position Rolotlve to the 
Detector at the Origin of the 
Sample Space Is (0,0,-a) 



FIGURE 12 



WO 99/26078 



PCT/US98/24Z20 



15/15 




cr 
o 
to 

CO 





UJ 








< 


■c 










1- 


<> 






UJ 



TKTT17DM A TTAXf A 1 Cn ntm/Mif 

IN IhJKJN A 1 1UNAL bEARCH REPORT 


Inter* ^vtl AppUcatlon No 
PfT/lK 0A/94??ft 


A. CLASSIFICATION OF SUBJECT MATTER 

IPC 6 G01R33/24 G01R33/16 




According to International Potont Classification (IPC) cx to both national classification and IPC 




B. FIELDS SEARCHED 


Mifwnum documentation eearched (classification system loflowoo Dy classification s/mDob) 

IPC 6 G01R 



Documentation searched other than rrHnimvm documentation to trie extent thai such documents aro included in me Ho? da searched 



ctecirofuc dala base consumed during the international search (name ol data base end. wfter© practical search terms used) 



C. DOCUMENTS C0NS1DEAEO TO 8€ RELEVANT 



Category " 


Citation ol document. «Mth indication, where appropriate, of the relevant passages 


Relevant ig c*airo No. 


X 


US 4 899 109 A (TROPP JAMES S ET AL) 


1-9, 




6 February 1990 


13-18, 






32,34, 






35,37, 






40-46 




see column 1, line 54 - column 1. line 68 






see column 4, line 27 - column 4, line 34 






see column 7, line 37 - column 9, line 36 






see column 13, line 18 - column 14, line 






21; figures 1,9,10 




A 




22,24. 






28-31 




-/-- 





0 



Fun her documents ore teJod in lt»o continual ion ol box C. 



0 



Patent lamry members &m lietod in annex. 



* Special categories ol cited documents : 

*A* document dehn:ng the general state ol the art which is not 

considered to be of particular *©far»anco 
*H" oerl** document Out pjD»3tiod on or «tter tha tate motional 

filing dato 

X* doc iment which rosy throw doubts on priority clalm<«) oc 
which is cited to establish the publication dale ol another 
elation or other special roaion (as specified) 

*0* documont raterring to an ord diacloeuro. use. exhibition or 
other means 

P* document putoiahod poor to Ihe rtemafonad Ming dote but 
rater than tho pnonty date cieimod 



*T" Utei document published ater rhe intomBDonal Ming, date 
or priority date and not fei conttct wtm the eppacatwn but 
cfiod to understand the principle or moon/ undertyina. tho 
invention 

'X~ document ol particular relevance: ihe claimed invention 
cannot too considered novel or cannot be considered to 
tnvorvo en inventive step when the document ra taken alono 

*Y* document ol particular letavanco: the cioffned invention 

cannot bo considered to iwvorve en inwjntrvo slop whoo the 
etoevmot* <c combined w*h one or more other auch docu- 
ments, euch combination being obvioue to a por aon cViftod 
in the art 

*&~ document member of the eat 



Dote ol ihe ectue) completion of the trtemetional anarch 

18 February 1999 


Oate ol mtilng of the international search report 

09/03/1999 


Name and making address of the tSA 

European Potent OfSce. P.O. 5818 Patent taan 2 
NL- 72901 IV Rijsw^lt 
Tot. (Oi ro)340>2040. T*. 31 651 epo rH. 
Fax: 4^31.70} 3«0^016 


Authorised officer 

Lersch, W 



Form PC T.1SA/2 1 0 {mcmm] ahMt) \±*r 1997) 



page 1 of 2 



INTERNATIONAL SEARCH REPORT 



litfori. .nel Application No 

PCT/US 98/24220 



C^Conilnu»t/on> DOCUMENTS CONSIDERED TO BE RELEVANT 



C«!eoo*y * 


^.iiuton o, ^xKrlui <uni fc wan kio *C-o u »*t t o u b py* ojjf iiiiij . ui n*tr »x*w«#im payyayi> 


fotoVAnt to cJftan No. 


Y 
A 


^ 1 t FT At • "Throo-Ri mane inna 1 Manriinn 
v. L J CI ML.. 1 II* t?tr u 1 Ultrllo I v»l la l l tajjp 1 uy 






UI vile OLd-LlC rldyiltLlL r ItrlU lit^lVJv uut 


17 lfi 




niitnan rieao 


1? 14 




magnfttt Rf^nNAwrr tn MFniriNF 


15 17 










too t ho uhnlo Hnrnmont 
^> tr" Lilt? vviiuic uuLUUicu l 




A 






A 


US S 073 858 A (MILLS RANDELL L) 


10-12, 




17 December 1991 


19.20, 




cited in tne application 














1ft 1Q 




see column j, line *4\j cuiujijii line tu 






see column o, line 00 column line 






see column iu, line 0 column /ii, ime d/ 






see column 26 > line 5 - column 27, line 






bl; figures Z-4,o,/ 




X 


KWIAT D ET AL: "A DECOUPLED COIL DETECTOR 


43 




ARRAY FOR FAST IMAGE ACQUISITION IN 






MAGNETIC RESONANCE IMAGING 






MEDICAL PHYSICS, 






vol. 18, no. 2, 1 March 1991, pages 






251-2ob, aP0002Zoo11 






see the whole document 




A 




1.17-19, 






21,32, 






11 -A 1 




md^ucc t 1 rr ai . * mat iurTi r 
nUorltK l 0 tl AL I nHliHLIIv. 


1 17-10 

1,1/17, 




CIKPFPTIRII T TY MrA^IIPFMFNT IKTNf A 






nniiRi F-nANTP tagging f nnr^ <;FfniFNf F" 


17-41 41 




MAGNETIC RESONANCE IN MEDICINE, 






vol. 18, no. 1, 1 March 1991, pages 






oct_occ Ypnnn?OQAii? 




A 


tic c i n A A ^TA&FUTr7 PfflTft M FT A1 ^ 


1 17-1Q 
1,1/ 1^ > 




l / nay i 


1? 






37-41,43 




see column 1, line 33 - column 1, line 52 






see column 3, line 15 - column 5, line 58; 






figures 1-3 




A 


P.T.CALLAGHAN: "Principles of Nuclear 


1.17-19, 




Magnetic Resonance Microscopy 0 


32, 




1995 , CLARENDON PRESS . OXFORD 


37-41,43 




XP002093385 






see page 208 - page 209 





Fwn PCI /ISA/210 (oantkntmiOA d frcCOftd W>*«f) puy 19V2) 



page 2 of 2 



INTERNATIONAL SEARCH REPORT 

..Kcx-mrlkon on patont !»fnfljf m«mb»r» 



In tot xv»1 Application Ho 

PCT/US 98/24220 



Patent document 




Publication 


Patent family 


Publication 


cited in search report 




date 


members} 


date 


US 4899109 


A 


06-02-1990 


EP 0434870 A 


03-07-1991 


US 5073858 


A 


17-12-1991 


US 4969469 A 


13-11-1990 


US 5313164 


A 


17-05-1994 


NONE 





Food PC MSA* IO {prtatfl t**nt f «nn««J IMy 1993) 



